From 51f1710d7eac046de2ba4cc580ceb030e37f57a7 Mon Sep 17 00:00:00 2001 From: Andreas Noack Date: Fri, 31 Jan 2020 08:16:33 +0100 Subject: [PATCH] Make SuiteSparse wrapper thread safe. (#34546) Fixes #34500 --- stdlib/SuiteSparse/src/cholmod.jl | 156 +++++++++++++++------------- stdlib/SuiteSparse/src/spqr.jl | 6 +- stdlib/SuiteSparse/test/cholmod.jl | 12 +-- stdlib/SuiteSparse/test/runtests.jl | 26 ++++- stdlib/SuiteSparse/test/threads.jl | 20 ++++ 5 files changed, 135 insertions(+), 85 deletions(-) create mode 100644 stdlib/SuiteSparse/test/threads.jl diff --git a/stdlib/SuiteSparse/src/cholmod.jl b/stdlib/SuiteSparse/src/cholmod.jl index 7b66ffcdfc32e..699de4398e70c 100644 --- a/stdlib/SuiteSparse/src/cholmod.jl +++ b/stdlib/SuiteSparse/src/cholmod.jl @@ -40,15 +40,15 @@ include("cholmod_h.jl") const CHOLMOD_MIN_VERSION = v"2.1.1" -const common_struct = Vector{UInt8}() +const common_struct = Vector{Vector{UInt8}}() -const common_supernodal = Ref{Ptr{Cint}}() -const common_final_ll = Ref{Ptr{Cint}}() -const common_print = Ref{Ptr{Cint}}() -const common_itype = Ref{Ptr{Cint}}() -const common_dtype = Ref{Ptr{Cint}}() -const common_nmethods = Ref{Ptr{Cint}}() -const common_postorder = Ref{Ptr{Cint}}() +const common_supernodal = Vector{Ptr{Cint}}() +const common_final_ll = Vector{Ptr{Cint}}() +const common_print = Vector{Ptr{Cint}}() +const common_itype = Vector{Ptr{Cint}}() +const common_dtype = Vector{Ptr{Cint}}() +const common_nmethods = Vector{Ptr{Cint}}() +const common_postorder = Vector{Ptr{Cint}}() ### These offsets are defined in SuiteSparse_wrapper.c const common_size = ccall((:jl_cholmod_common_size,:libsuitesparse_wrapper),Int,()) @@ -149,20 +149,30 @@ function __init__() ### Initiate CHOLMOD ### common_struct controls the type of factorization and keeps pointers - ### to temporary memory. - resize!(common_struct, common_size) - fill!(common_struct, 0xff) - - common_supernodal[] = pointer(common_struct, cholmod_com_offsets[4] + 1) - common_final_ll[] = pointer(common_struct, cholmod_com_offsets[7] + 1) - common_print[] = pointer(common_struct, cholmod_com_offsets[13] + 1) - common_itype[] = pointer(common_struct, cholmod_com_offsets[18] + 1) - common_dtype[] = pointer(common_struct, cholmod_com_offsets[19] + 1) - common_nmethods[] = pointer(common_struct, cholmod_com_offsets[15] + 1) - common_postorder[] = pointer(common_struct, cholmod_com_offsets[17] + 1) - - start(common_struct) # initializes CHOLMOD - set_print_level(common_struct, 0) # no printing from CHOLMOD by default + ### to temporary memory. We need to manage a copy for each thread. + nt = Threads.nthreads() + resize!(common_struct , nt) + resize!(common_supernodal, nt) + resize!(common_final_ll , nt) + resize!(common_print , nt) + resize!(common_itype , nt) + resize!(common_dtype , nt) + resize!(common_nmethods , nt) + resize!(common_postorder , nt) + for i in 1:nt + common_struct[i] = fill(0xff, common_size) + + common_supernodal[i] = pointer(common_struct[i], cholmod_com_offsets[4] + 1) + common_final_ll[i] = pointer(common_struct[i], cholmod_com_offsets[7] + 1) + common_print[i] = pointer(common_struct[i], cholmod_com_offsets[13] + 1) + common_itype[i] = pointer(common_struct[i], cholmod_com_offsets[18] + 1) + common_dtype[i] = pointer(common_struct[i], cholmod_com_offsets[19] + 1) + common_nmethods[i] = pointer(common_struct[i], cholmod_com_offsets[15] + 1) + common_postorder[i] = pointer(common_struct[i], cholmod_com_offsets[17] + 1) + + start(common_struct[i]) # initializes CHOLMOD + set_print_level(common_struct[i], 0) # no printing from CHOLMOD by default + end # Register gc tracked allocator if CHOLMOD is new enough if current_version >= v"3.0.0" @@ -179,7 +189,7 @@ function __init__() end function set_print_level(cm::Vector{UInt8}, lev::Integer) - unsafe_store!(common_print[], lev) + unsafe_store!(common_print[Threads.threadid()], lev) end #################### @@ -403,32 +413,32 @@ Factor(FC::FactorComponent) = Factor(FC.F) function allocate_dense(m::Integer, n::Integer, d::Integer, ::Type{Tv}) where {Tv<:VTypes} Dense(ccall((@cholmod_name("allocate_dense"), :libcholmod), Ptr{C_Dense{Tv}}, (Csize_t, Csize_t, Csize_t, Cint, Ptr{Cvoid}), - m, n, d, xtyp(Tv), common_struct)) + m, n, d, xtyp(Tv), common_struct[Threads.threadid()])) end function free!(p::Ptr{C_Dense{Tv}}) where {Tv<:VTypes} @isok ccall((@cholmod_name("free_dense"), :libcholmod), Cint, (Ref{Ptr{C_Dense{Tv}}}, Ptr{Cvoid}), - p, common_struct) + p, common_struct[Threads.threadid()]) end function zeros(m::Integer, n::Integer, ::Type{Tv}) where Tv<:VTypes Dense(ccall((@cholmod_name("zeros"), :libcholmod), Ptr{C_Dense{Tv}}, (Csize_t, Csize_t, Cint, Ptr{UInt8}), - m, n, xtyp(Tv), common_struct)) + m, n, xtyp(Tv), common_struct[Threads.threadid()])) end zeros(m::Integer, n::Integer) = zeros(m, n, Float64) function ones(m::Integer, n::Integer, ::Type{Tv}) where Tv<:VTypes Dense(ccall((@cholmod_name("ones"), :libcholmod), Ptr{C_Dense{Tv}}, (Csize_t, Csize_t, Cint, Ptr{UInt8}), - m, n, xtyp(Tv), common_struct)) + m, n, xtyp(Tv), common_struct[Threads.threadid()])) end ones(m::Integer, n::Integer) = ones(m, n, Float64) function eye(m::Integer, n::Integer, ::Type{Tv}) where Tv<:VTypes Dense(ccall((@cholmod_name("eye"), :libcholmod), Ptr{C_Dense{Tv}}, (Csize_t, Csize_t, Cint, Ptr{UInt8}), - m, n, xtyp(Tv), common_struct)) + m, n, xtyp(Tv), common_struct[Threads.threadid()])) end eye(m::Integer, n::Integer) = eye(m, n, Float64) eye(n::Integer) = eye(n, n, Float64) @@ -436,13 +446,13 @@ eye(n::Integer) = eye(n, n, Float64) function copy(A::Dense{Tv}) where Tv<:VTypes Dense(ccall((@cholmod_name("copy_dense"), :libcholmod), Ptr{C_Dense{Tv}}, (Ptr{C_Dense{Tv}}, Ptr{UInt8}), - A, common_struct)) + A, common_struct[Threads.threadid()])) end function sort!(S::Sparse{Tv}) where Tv<:VTypes @isok ccall((@cholmod_name("sort"), :libcholmod), Cint, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - S, common_struct) + S, common_struct[Threads.threadid()]) return S end @@ -458,14 +468,14 @@ function norm_dense(D::Dense{Tv}, p::Integer) where Tv<:VTypes end ccall((@cholmod_name("norm_dense"), :libcholmod), Cdouble, (Ptr{C_Dense{Tv}}, Cint, Ptr{UInt8}), - D, p, common_struct) + D, p, common_struct[Threads.threadid()]) end ### cholmod_check.h ### function check_dense(A::Dense{Tv}) where Tv<:VTypes ccall((@cholmod_name("check_dense"), :libcholmod), Cint, (Ptr{C_Dense{Tv}}, Ptr{UInt8}), - pointer(A), common_struct) != 0 + pointer(A), common_struct[Threads.threadid()]) != 0 end # Non-Dense wrappers @@ -477,40 +487,40 @@ function allocate_sparse(nrow::Integer, ncol::Integer, nzmax::Integer, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), nrow, ncol, nzmax, sorted, - packed, stype, xtyp(Tv), common_struct)) + packed, stype, xtyp(Tv), common_struct[Threads.threadid()])) end function free!(ptr::Ptr{C_Sparse{Tv}}) where Tv<:VTypes @isok ccall((@cholmod_name("free_sparse"), :libcholmod), Cint, (Ref{Ptr{C_Sparse{Tv}}}, Ptr{UInt8}), - ptr, common_struct) + ptr, common_struct[Threads.threadid()]) end function free!(ptr::Ptr{C_Factor{Tv}}) where Tv<:VTypes # Warning! Important that finalizer doesn't modify the global Common struct. @isok ccall((@cholmod_name("free_factor"), :libcholmod), Cint, (Ref{Ptr{C_Factor{Tv}}}, Ptr{Cvoid}), - ptr, common_struct) + ptr, common_struct[Threads.threadid()]) end function aat(A::Sparse{Tv}, fset::Vector{SuiteSparse_long}, mode::Integer) where Tv<:VRealTypes Sparse(ccall((@cholmod_name("aat"), :libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Ptr{SuiteSparse_long}, Csize_t, Cint, Ptr{UInt8}), - A, fset, length(fset), mode, common_struct)) + A, fset, length(fset), mode, common_struct[Threads.threadid()])) end function sparse_to_dense(A::Sparse{Tv}) where Tv<:VTypes Dense(ccall((@cholmod_name("sparse_to_dense"),:libcholmod), Ptr{C_Dense{Tv}}, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - A, common_struct)) + A, common_struct[Threads.threadid()])) end function dense_to_sparse(D::Dense{Tv}, ::Type{SuiteSparse_long}) where Tv<:VTypes Sparse(ccall((@cholmod_name("dense_to_sparse"),:libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Dense{Tv}}, Cint, Ptr{UInt8}), - D, true, common_struct)) + D, true, common_struct[Threads.threadid()])) end function factor_to_sparse!(F::Factor{Tv}) where Tv<:VTypes @@ -519,88 +529,88 @@ function factor_to_sparse!(F::Factor{Tv}) where Tv<:VTypes Sparse(ccall((@cholmod_name("factor_to_sparse"),:libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Factor{Tv}}, Ptr{UInt8}), - F, common_struct)) + F, common_struct[Threads.threadid()])) end function change_factor!(F::Factor{Tv}, to_ll::Bool, to_super::Bool, to_packed::Bool, to_monotonic::Bool) where Tv<:VTypes @isok ccall((@cholmod_name("change_factor"),:libcholmod), Cint, (Cint, Cint, Cint, Cint, Cint, Ptr{C_Factor{Tv}}, Ptr{UInt8}), - xtyp(Tv), to_ll, to_super, to_packed, to_monotonic, F, common_struct) + xtyp(Tv), to_ll, to_super, to_packed, to_monotonic, F, common_struct[Threads.threadid()]) end function check_sparse(A::Sparse{Tv}) where Tv<:VTypes ccall((@cholmod_name("check_sparse"),:libcholmod), Cint, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - A, common_struct) != 0 + A, common_struct[Threads.threadid()]) != 0 end function check_factor(F::Factor{Tv}) where Tv<:VTypes ccall((@cholmod_name("check_factor"),:libcholmod), Cint, (Ptr{C_Factor{Tv}}, Ptr{UInt8}), - F, common_struct) != 0 + F, common_struct[Threads.threadid()]) != 0 end function nnz(A::Sparse{Tv}) where Tv<:VTypes ccall((@cholmod_name("nnz"),:libcholmod), Int, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - A, common_struct) + A, common_struct[Threads.threadid()]) end function speye(m::Integer, n::Integer, ::Type{Tv}) where Tv<:VTypes Sparse(ccall((@cholmod_name("speye"), :libcholmod), Ptr{C_Sparse{Tv}}, (Csize_t, Csize_t, Cint, Ptr{UInt8}), - m, n, xtyp(Tv), common_struct)) + m, n, xtyp(Tv), common_struct[Threads.threadid()])) end function spzeros(m::Integer, n::Integer, nzmax::Integer, ::Type{Tv}) where Tv<:VTypes Sparse(ccall((@cholmod_name("spzeros"), :libcholmod), Ptr{C_Sparse{Tv}}, (Csize_t, Csize_t, Csize_t, Cint, Ptr{UInt8}), - m, n, nzmax, xtyp(Tv), common_struct)) + m, n, nzmax, xtyp(Tv), common_struct[Threads.threadid()])) end function transpose_(A::Sparse{Tv}, values::Integer) where Tv<:VTypes Sparse(ccall((@cholmod_name("transpose"),:libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Cint, Ptr{UInt8}), - A, values, common_struct)) + A, values, common_struct[Threads.threadid()])) end function copy(F::Factor{Tv}) where Tv<:VTypes Factor(ccall((@cholmod_name("copy_factor"),:libcholmod), Ptr{C_Factor{Tv}}, (Ptr{C_Factor{Tv}}, Ptr{UInt8}), - F, common_struct)) + F, common_struct[Threads.threadid()])) end function copy(A::Sparse{Tv}) where Tv<:VTypes Sparse(ccall((@cholmod_name("copy_sparse"),:libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - A, common_struct)) + A, common_struct[Threads.threadid()])) end function copy(A::Sparse{Tv}, stype::Integer, mode::Integer) where Tv<:VRealTypes Sparse(ccall((@cholmod_name("copy"),:libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Cint, Cint, Ptr{UInt8}), - A, stype, mode, common_struct)) + A, stype, mode, common_struct[Threads.threadid()])) end ### cholmod_check.h ### function print_sparse(A::Sparse{Tv}, name::String) where Tv<:VTypes isascii(name) || error("non-ASCII name: $name") - set_print_level(common_struct, 3) + set_print_level(common_struct[Threads.threadid()], 3) @isok ccall((@cholmod_name("print_sparse"),:libcholmod), Cint, (Ptr{C_Sparse{Tv}}, Ptr{UInt8}, Ptr{UInt8}), - A, name, common_struct) + A, name, common_struct[Threads.threadid()]) nothing end function print_factor(F::Factor{Tv}, name::String) where Tv<:VTypes - set_print_level(common_struct, 3) + set_print_level(common_struct[Threads.threadid()], 3) @isok ccall((@cholmod_name("print_factor"),:libcholmod), Cint, (Ptr{C_Factor{Tv}}, Ptr{UInt8}, Ptr{UInt8}), - F, name, common_struct) + F, name, common_struct[Threads.threadid()]) nothing end @@ -617,7 +627,7 @@ function ssmult(A::Sparse{Tv}, B::Sparse{Tv}, stype::Integer, (Ptr{C_Sparse{Tv}}, Ptr{C_Sparse{Tv}}, Cint, Cint, Cint, Ptr{UInt8}), A, B, stype, values, - sorted, common_struct)) + sorted, common_struct[Threads.threadid()])) end function norm_sparse(A::Sparse{Tv}, norm::Integer) where Tv<:VTypes @@ -626,14 +636,14 @@ function norm_sparse(A::Sparse{Tv}, norm::Integer) where Tv<:VTypes end ccall((@cholmod_name("norm_sparse"), :libcholmod), Cdouble, (Ptr{C_Sparse{Tv}}, Cint, Ptr{UInt8}), - A, norm, common_struct) + A, norm, common_struct[Threads.threadid()]) end function horzcat(A::Sparse{Tv}, B::Sparse{Tv}, values::Bool) where Tv<:VRealTypes Sparse(ccall((@cholmod_name("horzcat"), :libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Ptr{C_Sparse{Tv}}, Cint, Ptr{UInt8}), - A, B, values, common_struct)) + A, B, values, common_struct[Threads.threadid()])) end function scale!(S::Dense{Tv}, scale::Integer, A::Sparse{Tv}) where Tv<:VRealTypes @@ -662,7 +672,7 @@ function scale!(S::Dense{Tv}, scale::Integer, A::Sparse{Tv}) where Tv<:VRealType sA = unsafe_load(pointer(A)) @isok ccall((@cholmod_name("scale"),:libcholmod), Cint, (Ptr{C_Dense{Tv}}, Cint, Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - S, scale, A, common_struct) + S, scale, A, common_struct[Threads.threadid()]) A end @@ -678,7 +688,7 @@ function sdmult!(A::Sparse{Tv}, transpose::Bool, (Ptr{C_Sparse{Tv}}, Cint, Ref{ComplexF64}, Ref{ComplexF64}, Ptr{C_Dense{Tv}}, Ptr{C_Dense{Tv}}, Ptr{UInt8}), - A, transpose, α, β, X, Y, common_struct) + A, transpose, α, β, X, Y, common_struct[Threads.threadid()]) Y end @@ -686,7 +696,7 @@ function vertcat(A::Sparse{Tv}, B::Sparse{Tv}, values::Bool) where Tv<:VRealType Sparse(ccall((@cholmod_name("vertcat"), :libcholmod), Ptr{C_Sparse{Tv}}, (Ptr{C_Sparse{Tv}}, Ptr{C_Sparse{Tv}}, Cint, Ptr{UInt8}), - A, B, values, common_struct)) + A, B, values, common_struct[Threads.threadid()])) end function symmetry(A::Sparse{Tv}, option::Integer) where Tv<:VTypes @@ -698,7 +708,7 @@ function symmetry(A::Sparse{Tv}, option::Integer) where Tv<:VTypes (Ptr{C_Sparse{Tv}}, Cint, Ptr{SuiteSparse_long}, Ptr{SuiteSparse_long}, Ptr{SuiteSparse_long}, Ptr{SuiteSparse_long}, Ptr{UInt8}), A, option, xmatched, pmatched, - nzoffdiag, nzdiag, common_struct) + nzoffdiag, nzdiag, common_struct[Threads.threadid()]) rv, xmatched[], pmatched[], nzoffdiag[], nzdiag[] end @@ -751,7 +761,7 @@ function solve(sys::Integer, F::Factor{Tv}, B::Dense{Tv}) where Tv<:VTypes end Dense(ccall((@cholmod_name("solve"),:libcholmod), Ptr{C_Dense{Tv}}, (Cint, Ptr{C_Factor{Tv}}, Ptr{C_Dense{Tv}}, Ptr{UInt8}), - sys, F, B, common_struct)) + sys, F, B, common_struct[Threads.threadid()])) end function spsolve(sys::Integer, F::Factor{Tv}, B::Sparse{Tv}) where Tv<:VTypes @@ -762,7 +772,7 @@ function spsolve(sys::Integer, F::Factor{Tv}, B::Sparse{Tv}) where Tv<:VTypes Sparse(ccall((@cholmod_name("spsolve"),:libcholmod), Ptr{C_Sparse{Tv}}, (Cint, Ptr{C_Factor{Tv}}, Ptr{C_Sparse{Tv}}, Ptr{UInt8}), - sys, F, B, common_struct)) + sys, F, B, common_struct[Threads.threadid()])) end # Autodetects the types @@ -770,7 +780,7 @@ function read_sparse(file::Libc.FILE, ::Type{SuiteSparse_long}) ptr = ccall((@cholmod_name("read_sparse"), :libcholmod), Ptr{C_Sparse{Cvoid}}, (Ptr{Cvoid}, Ptr{UInt8}), - file.ptr, common_struct) + file.ptr, common_struct[Threads.threadid()]) if ptr == C_NULL throw(ArgumentError("sparse matrix construction failed. Check that input file is valid.")) end @@ -1272,14 +1282,14 @@ function fact_(A::Sparse{<:VTypes}, cm::Array{UInt8}; sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) if !postorder - unsafe_store!(common_postorder[], 0) + unsafe_store!(common_postorder[Threads.threadid()], 0) end if perm === nothing || isempty(perm) # TODO: deprecate empty perm F = analyze(A, cm) else # user permutation provided if userperm_only # use perm even if it is worse than AMD - unsafe_store!(common_nmethods[], 1) + unsafe_store!(common_nmethods[Threads.threadid()], 1) end F = analyze_p(A, SuiteSparse_long[p-1 for p in perm], cm) end @@ -1290,10 +1300,10 @@ end function cholesky!(F::Factor{Tv}, A::Sparse{Tv}; shift::Real=0.0, check::Bool = true) where Tv # Makes it an LLt - unsafe_store!(common_final_ll[], 1) + unsafe_store!(common_final_ll[Threads.threadid()], 1) # Compute the numerical factorization - factorize_p!(A, shift, F, common_struct) + factorize_p!(A, shift, F, common_struct[Threads.threadid()]) check && (issuccess(F) || throw(LinearAlgebra.PosDefException(1))) return F @@ -1326,7 +1336,7 @@ cholesky!(F::Factor, A::Union{SparseMatrixCSC{T}, function cholesky(A::Sparse; shift::Real=0.0, check::Bool = true, perm::Union{Nothing,AbstractVector{SuiteSparse_long}}=nothing) - cm = defaults(common_struct) + cm = defaults(common_struct[Threads.threadid()]) set_print_level(cm, 0) # Compute the symbolic factorization @@ -1450,7 +1460,7 @@ cholesky(A::Union{SparseMatrixCSC{T}, SparseMatrixCSC{Complex{T}}, function ldlt!(F::Factor{Tv}, A::Sparse{Tv}; shift::Real=0.0, check::Bool = true) where Tv - cm = defaults(common_struct) + cm = defaults(common_struct[Threads.threadid()]) set_print_level(cm, 0) # Makes it an LDLt @@ -1490,13 +1500,13 @@ ldlt!(F::Factor, A::Union{SparseMatrixCSC{T}, function ldlt(A::Sparse; shift::Real=0.0, check::Bool = true, perm::Union{Nothing,AbstractVector{SuiteSparse_long}}=nothing) - cm = defaults(common_struct) + cm = defaults(common_struct[Threads.threadid()]) set_print_level(cm, 0) # Makes it an LDLt - unsafe_store!(common_final_ll[], 0) + unsafe_store!(common_final_ll[Threads.threadid()], 0) # Really make sure it's an LDLt by avoiding supernodal factorization - unsafe_store!(common_supernodal[], 0) + unsafe_store!(common_supernodal[Threads.threadid()], 0) # Compute the symbolic factorization F = fact_(A, cm; perm = perm) @@ -1571,7 +1581,7 @@ function lowrankupdowndate!(F::Factor{Tv}, C::Sparse{Tv}, update::Cint) where Tv end @isok ccall((@cholmod_name("updown"), :libcholmod), Cint, (Cint, Ptr{C_Sparse{Tv}}, Ptr{C_Factor{Tv}}, Ptr{Cvoid}), - update, C, F, common_struct) + update, C, F, common_struct[Threads.threadid()]) F end diff --git a/stdlib/SuiteSparse/src/spqr.jl b/stdlib/SuiteSparse/src/spqr.jl index d5bc42562a8c7..983f7ba37e017 100644 --- a/stdlib/SuiteSparse/src/spqr.jl +++ b/stdlib/SuiteSparse/src/spqr.jl @@ -63,7 +63,7 @@ function _qr!(ordering::Integer, tol::Real, econ::Integer, getCTX::Integer, H, # m-by-nh Householder vectors HPinv, # size m row permutation HTau, # 1-by-nh Householder coefficients - CHOLMOD.common_struct) # /* workspace and parameters */ + CHOLMOD.common_struct[Threads.threadid()]) # /* workspace and parameters */ if rnk < 0 error("Sparse QR factorization failed") @@ -82,7 +82,7 @@ function _qr!(ordering::Integer, tol::Real, econ::Integer, getCTX::Integer, # the common struct is updated ccall((:cholmod_l_free, :libcholmod), Cvoid, (Csize_t, Cint, Ptr{CHOLMOD.SuiteSparse_long}, Ptr{Cvoid}), - n, sizeof(CHOLMOD.SuiteSparse_long), e, CHOLMOD.common_struct) + n, sizeof(CHOLMOD.SuiteSparse_long), e, CHOLMOD.common_struct[Threads.threadid()]) end hpinv = HPinv[] if hpinv == C_NULL @@ -97,7 +97,7 @@ function _qr!(ordering::Integer, tol::Real, econ::Integer, getCTX::Integer, # the common struct is updated ccall((:cholmod_l_free, :libcholmod), Cvoid, (Csize_t, Cint, Ptr{CHOLMOD.SuiteSparse_long}, Ptr{Cvoid}), - m, sizeof(CHOLMOD.SuiteSparse_long), hpinv, CHOLMOD.common_struct) + m, sizeof(CHOLMOD.SuiteSparse_long), hpinv, CHOLMOD.common_struct[Threads.threadid()]) end return rnk, _E, _HPinv diff --git a/stdlib/SuiteSparse/test/cholmod.jl b/stdlib/SuiteSparse/test/cholmod.jl index 5c3be022daa1a..c19c3726dd979 100644 --- a/stdlib/SuiteSparse/test/cholmod.jl +++ b/stdlib/SuiteSparse/test/cholmod.jl @@ -234,7 +234,7 @@ end @testset "illegal dtype (for now but should be supported at some point)" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Cvoid}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) puint = convert(Ptr{UInt32}, p) unsafe_store!(puint, CHOLMOD.SINGLE, 3*div(sizeof(Csize_t), 4) + 5*div(sizeof(Ptr{Cvoid}), 4) + 4) @test_throws CHOLMOD.CHOLMODException CHOLMOD.Sparse(p) @@ -243,7 +243,7 @@ end @testset "illegal dtype" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Cvoid}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) puint = convert(Ptr{UInt32}, p) unsafe_store!(puint, 5, 3*div(sizeof(Csize_t), 4) + 5*div(sizeof(Ptr{Cvoid}), 4) + 4) @test_throws CHOLMOD.CHOLMODException CHOLMOD.Sparse(p) @@ -252,7 +252,7 @@ end @testset "illegal xtype" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Cvoid}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) puint = convert(Ptr{UInt32}, p) unsafe_store!(puint, 3, 3*div(sizeof(Csize_t), 4) + 5*div(sizeof(Ptr{Cvoid}), 4) + 3) @test_throws CHOLMOD.CHOLMODException CHOLMOD.Sparse(p) @@ -261,7 +261,7 @@ end @testset "illegal itype I" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Cvoid}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) puint = convert(Ptr{UInt32}, p) unsafe_store!(puint, CHOLMOD.INTLONG, 3*div(sizeof(Csize_t), 4) + 5*div(sizeof(Ptr{Cvoid}), 4) + 2) @test_throws CHOLMOD.CHOLMODException CHOLMOD.Sparse(p) @@ -270,7 +270,7 @@ end @testset "illegal itype II" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Cvoid}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) puint = convert(Ptr{UInt32}, p) unsafe_store!(puint, 5, 3*div(sizeof(Csize_t), 4) + 5*div(sizeof(Ptr{Cvoid}), 4) + 2) @test_throws CHOLMOD.CHOLMODException CHOLMOD.Sparse(p) @@ -321,7 +321,7 @@ end @testset "test free!" begin p = ccall((:cholmod_l_allocate_sparse, :libcholmod), Ptr{CHOLMOD.C_Sparse{Float64}}, (Csize_t, Csize_t, Csize_t, Cint, Cint, Cint, Cint, Ptr{Cvoid}), - 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct) + 1, 1, 1, true, true, 0, CHOLMOD.REAL, CHOLMOD.common_struct[Threads.threadid()]) @test CHOLMOD.free!(p) end diff --git a/stdlib/SuiteSparse/test/runtests.jl b/stdlib/SuiteSparse/test/runtests.jl index d317023c6c84d..cde54e9488818 100644 --- a/stdlib/SuiteSparse/test/runtests.jl +++ b/stdlib/SuiteSparse/test/runtests.jl @@ -4,7 +4,27 @@ using Test, Random using SuiteSparse, LinearAlgebra, SparseArrays if Base.USE_GPL_LIBS - include("umfpack.jl") - include("cholmod.jl") - include("spqr.jl") + include("umfpack.jl") + include("cholmod.jl") + include("spqr.jl") + + # Test multithreaded execution + let p, cmd = `$(Base.julia_cmd()) --depwarn=error --startup-file=no threads.jl` + # test both nthreads==1 and nthreads>1. spawn a process to test whichever + # case we are not running currently. + other_nthreads = Threads.nthreads() == 1 ? 4 : 1 + p = run( + pipeline( + setenv( + cmd, + "JULIA_NUM_THREADS" => other_nthreads, + dir=@__DIR__()), + stdout = stdout, + stderr = stderr), + wait = false) + include("threads.jl") + if !success(p) + error("SuiteSparse threads test failed with nthreads == $other_nthreads") + end + end end diff --git a/stdlib/SuiteSparse/test/threads.jl b/stdlib/SuiteSparse/test/threads.jl new file mode 100644 index 0000000000000..76b4acd540779 --- /dev/null +++ b/stdlib/SuiteSparse/test/threads.jl @@ -0,0 +1,20 @@ +using Test, LinearAlgebra, SparseArrays + +@testset "threaded SuiteSparse tests" begin + A = sprandn(200, 200, 0.2) + b = rand(200) + + function test(n::Integer) + _A = A[1:n, 1:n] + _b = b[1:n] + x = qr(_A) \ _b + return norm(x) + end + + res_threads = zeros(100) + Threads.@threads for i in 1:100 + res_threads[i] = test(i + 100) + end + + @test res_threads ≈ [test(i + 100) for i in 1:100] +end