From 5d3724a59255d3bdacd1eae146c20da317e30ccd Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Mon, 21 Oct 2024 17:05:38 +0200 Subject: [PATCH 1/8] Inline sparse-times-dense in-place multiplication --- src/linalg.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/linalg.jl b/src/linalg.jl index 2251c3ae..c76f26a6 100644 --- a/src/linalg.jl +++ b/src/linalg.jl @@ -28,11 +28,11 @@ for op ∈ (:+, :-) end end -LinearAlgebra.generic_matmatmul!(C::StridedMatrix, tA, tB, A::SparseMatrixCSCUnion, B::DenseMatrixUnion, _add::MulAddMul) = +@inline LinearAlgebra.generic_matmatmul!(C::StridedMatrix, tA, tB, A::SparseMatrixCSCUnion, B::DenseMatrixUnion, _add::MulAddMul) = spdensemul!(C, tA, tB, A, B, _add) -LinearAlgebra.generic_matmatmul!(C::StridedMatrix, tA, tB, A::SparseMatrixCSCUnion, B::AbstractTriangular, _add::MulAddMul) = +@inline LinearAlgebra.generic_matmatmul!(C::StridedMatrix, tA, tB, A::SparseMatrixCSCUnion, B::AbstractTriangular, _add::MulAddMul) = spdensemul!(C, tA, tB, A, B, _add) -LinearAlgebra.generic_matvecmul!(C::StridedVecOrMat, tA, A::SparseMatrixCSCUnion, B::DenseInputVector, _add::MulAddMul) = +@inline LinearAlgebra.generic_matvecmul!(C::StridedVecOrMat, tA, A::SparseMatrixCSCUnion, B::DenseInputVector, _add::MulAddMul) = spdensemul!(C, tA, 'N', A, B, _add) Base.@constprop :aggressive function spdensemul!(C, tA, tB, A, B, _add) From 91b0aa56a64136c920799533ab6a743d008118bf Mon Sep 17 00:00:00 2001 From: CyHan Date: Sat, 26 Oct 2024 22:29:08 +0800 Subject: [PATCH 2/8] doc: move solvers doc to `src\solvers.md` (#576) --- docs/src/index.md | 5 +++++ docs/src/solvers.md | 23 +++++++++++++++-------- 2 files changed, 20 insertions(+), 8 deletions(-) diff --git a/docs/src/index.md b/docs/src/index.md index 5ddae840..ceca6fb9 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -202,6 +202,11 @@ section of the standard library reference. | [`sprandn(m,n,d)`](@ref) | [`randn(m,n)`](@ref) | Creates a *m*-by-*n* random matrix (of density *d*) with iid non-zero elements distributed according to the standard normal (Gaussian) distribution. | | [`sprandn(rng,m,n,d)`](@ref) | [`randn(rng,m,n)`](@ref) | Creates a *m*-by-*n* random matrix (of density *d*) with iid non-zero elements generated with the `rng` random number generator | + +```@meta +DocTestSetup = nothing +``` + # [SparseArrays API](@id stdlib-sparse-arrays) ```@docs diff --git a/docs/src/solvers.md b/docs/src/solvers.md index e633be9d..4763ffcb 100644 --- a/docs/src/solvers.md +++ b/docs/src/solvers.md @@ -4,7 +4,16 @@ DocTestSetup = :(using LinearAlgebra, SparseArrays) ``` -Sparse matrix solvers call functions from [SuiteSparse](http://suitesparse.com). The following factorizations are available: +## [Sparse Linear Algebra](@id stdlib-sparse-linalg) + +Sparse matrix solvers call functions from [SuiteSparse](http://suitesparse.com). + +The following factorizations are available: + +1. [`cholesky`](@ref SparseArrays.CHOLMOD.cholesky) +2. [`ldlt`](@ref SparseArrays.CHOLMOD.ldlt) +3. [`lu`](@ref SparseArrays.UMFPACK.lu) +4. [`qr`](@ref SparseArrays.SPQR.qr) | Type | Description | |:--------------------------------- |:--------------------------------------------- | @@ -12,18 +21,16 @@ Sparse matrix solvers call functions from [SuiteSparse](http://suitesparse.com). | `UMFPACK.UmfpackLU` | LU factorization | | `SPQR.QRSparse` | QR factorization | -Other solvers such as [Pardiso.jl](https://github.com/JuliaSparse/Pardiso.jl/) are available as external packages. [Arpack.jl](https://julialinearalgebra.github.io/Arpack.jl/stable/) provides `eigs` and `svds` for iterative solution of eigensystems and singular value decompositions. +Other solvers such as [Pardiso.jl](https://github.com/JuliaSparse/Pardiso.jl/) are available +as external packages. [Arpack.jl](https://julialinearalgebra.github.io/Arpack.jl/stable/) +provides `eigs` and `svds` for iterative solution of eigensystems and singular value +decompositions. These factorizations are described in more detail in the [`Linear Algebra`](https://docs.julialang.org/en/v1/stdlib/LinearAlgebra/) section of the manual: -1. [`cholesky`](@ref SparseArrays.CHOLMOD.cholesky) -2. [`ldlt`](@ref SparseArrays.CHOLMOD.ldlt) -3. [`lu`](@ref SparseArrays.UMFPACK.lu) -4. [`qr`](@ref SparseArrays.SPQR.qr) - -```@docs +```@docs; canonical=false SparseArrays.CHOLMOD.cholesky SparseArrays.CHOLMOD.cholesky! SparseArrays.CHOLMOD.ldlt From 30fbfc6b0f50a9f1932248e0543ba4714c4e6f35 Mon Sep 17 00:00:00 2001 From: Dilum Aluthge Date: Wed, 23 Aug 2023 18:28:46 -0400 Subject: [PATCH 3/8] Test suite: activate a temp project if we need to install Aqua.jl during the test suite (#425) --- test/ambiguous.jl | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/test/ambiguous.jl b/test/ambiguous.jl index 4b81ae2f..c1d4ec58 100644 --- a/test/ambiguous.jl +++ b/test/ambiguous.jl @@ -5,6 +5,7 @@ original_depot_path = copy(Base.DEPOT_PATH) original_load_path = copy(Base.LOAD_PATH) original_env = copy(ENV) +original_project = Base.active_project() ### import Pkg @@ -13,6 +14,7 @@ import Pkg if Base.find_package("Aqua") === nothing @debug "Installing Aqua.jl for SparseArrays.jl tests" iob = IOBuffer() + Pkg.activate(; temp = true) try # TODO: make this version tie to compat in Project.toml # or do this another safer way @@ -103,4 +105,6 @@ end for (k, v) in pairs(original_env) ENV[k] = v end + +Base.set_active_project(original_project) ### From 546be180bbc79b99d0faba99864a19ecd1a031ab Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Sat, 26 Aug 2023 12:04:00 -0400 Subject: [PATCH 4/8] Fix docs conflict when building as part of full Julia docs (#430) Fix docs conflict --- docs/src/index.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/src/index.md b/docs/src/index.md index ceca6fb9..dc418088 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -7,7 +7,7 @@ DocTestSetup = :(using SparseArrays, LinearAlgebra) Julia has support for sparse vectors and [sparse matrices](https://en.wikipedia.org/wiki/Sparse_matrix) in the `SparseArrays` stdlib module. Sparse arrays are arrays that contain enough zeros that storing them in a special data structure leads to savings in space and execution time, compared to dense arrays. -External packages which implement different sparse storage types, multidimensional sparse arrays, and more can be found in [Noteworthy external packages](@ref man-csc) +External packages which implement different sparse storage types, multidimensional sparse arrays, and more can be found in [Noteworthy External Sparse Packages](@ref) ## [Compressed Sparse Column (CSC) Sparse Matrix Storage](@id man-csc) @@ -246,7 +246,7 @@ SparseArrays.ftranspose! ```@meta DocTestSetup = nothing ``` -# Noteworthy external packages +# Noteworthy External Sparse Packages Several other Julia packages provide sparse matrix implementations that should be mentioned: From 9b8cd143b0512239f6c3b513a9babe245985ff58 Mon Sep 17 00:00:00 2001 From: Jishnu Bhattacharya Date: Thu, 4 Apr 2024 08:24:39 +0530 Subject: [PATCH 5/8] SparseMatrixCSC constructor with a Tuple of Integers (#523) * SparseMatrixCSC constructor with a Tuple of Integers * SparseVector constructors --------- Co-authored-by: Viral B. Shah --- src/sparsematrix.jl | 2 ++ src/sparsevector.jl | 1 + test/sparsematrix_constructors_indexing.jl | 11 +++++++---- test/sparsevector.jl | 12 ++++++++---- 4 files changed, 18 insertions(+), 8 deletions(-) diff --git a/src/sparsematrix.jl b/src/sparsematrix.jl index 1703c870..4abb87ab 100644 --- a/src/sparsematrix.jl +++ b/src/sparsematrix.jl @@ -49,10 +49,12 @@ SparseMatrixCSC(m, n, colptr::ReadOnly, rowval::ReadOnly, nzval::Vector) = """ SparseMatrixCSC{Tv,Ti}(::UndefInitializer, m::Integer, n::Integer) + SparseMatrixCSC{Tv,Ti}(::UndefInitializer, (m,n)::NTuple{2,Integer}) Creates an empty sparse matrix with element type `Tv` and integer type `Ti` of size `m × n`. """ SparseMatrixCSC{Tv,Ti}(::UndefInitializer, m::Integer, n::Integer) where {Tv, Ti} = spzeros(Tv, Ti, m, n) +SparseMatrixCSC{Tv,Ti}(::UndefInitializer, mn::NTuple{2,Integer}) where {Tv, Ti} = spzeros(Tv, Ti, mn...) """ `FixedSparseCSC{Tv,Ti<:Integer} <: AbstractSparseMatrixCSC{Tv,Ti}` diff --git a/src/sparsevector.jl b/src/sparsevector.jl index 37b365ab..1cfe3b4f 100644 --- a/src/sparsevector.jl +++ b/src/sparsevector.jl @@ -50,6 +50,7 @@ SparseVector(n::Integer, nzind::Vector{Ti}, nzval::Vector{Tv}) where {Tv,Ti} = SparseVector{Tv,Ti}(n, nzind, nzval) SparseVector{Tv, Ti}(::UndefInitializer, n::Integer) where {Tv, Ti} = SparseVector{Tv, Ti}(n, Ti[], Tv[]) +SparseVector{Tv, Ti}(::UndefInitializer, (n,)::Tuple{Integer}) where {Tv, Ti} = SparseVector{Tv, Ti}(n, Ti[], Tv[]) """ `FixedSparseVector{Tv,Ti<:Integer} <: AbstractCompressedVector{Tv,Ti}` diff --git a/test/sparsematrix_constructors_indexing.jl b/test/sparsematrix_constructors_indexing.jl index 5aad5628..a6f1a75c 100644 --- a/test/sparsematrix_constructors_indexing.jl +++ b/test/sparsematrix_constructors_indexing.jl @@ -57,10 +57,13 @@ end @test sparse([1, 1, 2, 2, 2], [1, 2, 1, 2, 2], -1.0, 2, 2, *) == sparse([1, 1, 2, 2], [1, 2, 1, 2], [-1.0, -1.0, -1.0, 1.0], 2, 2) @test sparse(sparse(Int32.(1:5), Int32.(1:5), trues(5))') isa SparseMatrixCSC{Bool,Int32} # undef initializer - m = SparseMatrixCSC{Float32, Int16}(undef, 3, 4) - @test size(m) == (3, 4) - @test eltype(m) === Float32 - @test m == spzeros(3, 4) + sz = (3, 4) + for m in (SparseMatrixCSC{Float32, Int16}(undef, sz...), SparseMatrixCSC{Float32, Int16}(undef, sz), + similar(SparseMatrixCSC{Float32, Int16}, sz)) + @test size(m) == sz + @test eltype(m) === Float32 + @test m == spzeros(sz...) + end end @testset "spzeros for pattern creation (structural zeros)" begin diff --git a/test/sparsevector.jl b/test/sparsevector.jl index 44659939..b93a4963 100644 --- a/test/sparsevector.jl +++ b/test/sparsevector.jl @@ -219,10 +219,14 @@ end end @testset "Undef initializer" begin - v = SparseVector{Float32, Int16}(undef, 4) - @test size(v) == (4, ) - @test eltype(v) === Float32 - @test v == spzeros(Float32, 4) + sz = (4,) + for v in (SparseVector{Float32, Int16}(undef, sz), + SparseVector{Float32, Int16}(undef, sz...), + similar(SparseVector{Float32, Int16}, sz)) + @test size(v) == sz + @test eltype(v) === Float32 + @test v == spzeros(Float32, sz...) + end end end ### Element access From d2a80a6741d721a3dfebf42d5a12531133fb23c7 Mon Sep 17 00:00:00 2001 From: Ben Corbett <32752943+corbett5@users.noreply.github.com> Date: Thu, 22 Aug 2024 11:40:53 -0600 Subject: [PATCH 6/8] Change default QR tolerance to match SPQR (#557) SPQR uses just the double precision epsilon even for Float32. https://github.com/DrTimothyAldenDavis/SuiteSparse/blob/131471310ef0600b231b8fa7c10a55c3f70afbd9/SPQR/Source/spqr_tol.cpp#L29C6-L30C57 --- src/solvers/spqr.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/solvers/spqr.jl b/src/solvers/spqr.jl index d1a09f90..1b66f14b 100644 --- a/src/solvers/spqr.jl +++ b/src/solvers/spqr.jl @@ -146,7 +146,7 @@ Matrix{T}(Q::QRSparseQ) where {T} = lmul!(Q, Matrix{T}(I, size(Q, 1), min(size(Q # From SPQR manual p. 6 _default_tol(A::AbstractSparseMatrixCSC) = - 20*sum(size(A))*eps(real(eltype(A)))*maximum(norm(view(A, :, i)) for i in 1:size(A, 2)) + 20*sum(size(A))*eps()*maximum(norm(view(A, :, i)) for i in 1:size(A, 2)) """ qr(A::SparseMatrixCSC; tol=_default_tol(A), ordering=ORDERING_DEFAULT) -> QRSparse From 2d762b3593499ccf35e7ef0cc4baaaf3cdb20ce6 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" <744411+ViralBShah@users.noreply.github.com> Date: Sat, 2 Nov 2024 12:09:45 -0400 Subject: [PATCH 7/8] Manual commit for PR #550 to backport to 1.10 (#577) * Disable nested dissection (from #550) * Update CI for Julia 1.10 --------- Co-authored-by: Viral B. Shah --- .github/workflows/ci.yml | 23 +++++++++++++---------- src/solvers/cholmod.jl | 28 ++++++++++++++++++---------- 2 files changed, 31 insertions(+), 20 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index faeaf424..5a9fad11 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -22,24 +22,26 @@ jobs: fail-fast: false matrix: version: - - '~1.10.0-0' + - '1.10' os: - ubuntu-latest - - macOS-latest - windows-latest arch: - x64 - - x86 - exclude: + include: - os: macOS-latest + arch: aarch64 + version: '1.10' + - os: ubuntu-latest arch: x86 + version: '1.10' steps: - - uses: actions/checkout@v2 - - uses: julia-actions/setup-julia@v1 + - uses: actions/checkout@v4 + - uses: julia-actions/setup-julia@v2 with: version: ${{ matrix.version }} arch: ${{ matrix.arch }} - - uses: actions/cache@v1 + - uses: actions/cache@v4 env: cache-name: cache-artifacts with: @@ -55,17 +57,18 @@ jobs: env: SPARSEARRAYS_AQUA_TEST: true - uses: julia-actions/julia-processcoverage@v1 - - uses: codecov/codecov-action@v1 + - uses: codecov/codecov-action@v4 with: file: lcov.info + token: ${{ secrets.CODECOV_TOKEN }} + docs: runs-on: ubuntu-latest steps: - uses: actions/checkout@v2 - uses: julia-actions/setup-julia@latest with: - # version: '1.6' - version: 'nightly' + version: '1.10' - name: Generate docs run: | julia --color=yes -e 'write("Project.toml", replace(read("Project.toml", String), r"uuid = .*?\n" =>"uuid = \"3f01184e-e22b-5df5-ae63-d93ebab69eaf\"\n"));' diff --git a/src/solvers/cholmod.jl b/src/solvers/cholmod.jl index 731eaa3f..97ca426e 100644 --- a/src/solvers/cholmod.jl +++ b/src/solvers/cholmod.jl @@ -1258,25 +1258,33 @@ end ## Compute that symbolic factorization only function symbolic(A::Sparse{<:VTypes, Ti}; - perm::Union{Nothing,AbstractVector{<:Integer}}=nothing, - postorder::Bool=isnothing(perm)||isempty(perm), userperm_only::Bool=true) where Ti + perm::Union{Nothing,AbstractVector{<:Integer}}=nothing, + postorder::Bool=isnothing(perm)||isempty(perm), + userperm_only::Bool=true, + nested_dissection::Bool=false) where Ti sA = unsafe_load(pointer(A)) sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) - @cholmod_param postorder = postorder begin - if perm === nothing || isempty(perm) # TODO: deprecate empty perm - return analyze(A) - else # user permutation provided - if userperm_only # use perm even if it is worse than AMD - @cholmod_param nmethods = 1 begin + # The default is to just use AMD. Use nested dissection only if explicitly asked for. + # https://github.com/JuliaSparse/SparseArrays.jl/issues/548 + # https://github.com/DrTimothyAldenDavis/SuiteSparse/blob/26ababc7f3af725c5fb9168a1b94850eab74b666/CHOLMOD/Include/cholmod.h#L555-L574 + @cholmod_param nmethods = (nested_dissection ? 0 : 2) begin + @cholmod_param postorder = postorder begin + if perm === nothing || isempty(perm) # TODO: deprecate empty perm + return analyze(A) + else # user permutation provided + if userperm_only # use perm even if it is worse than AMD + @cholmod_param nmethods = 1 begin + return analyze_p(A, Ti[p-1 for p in perm]) + end + else return analyze_p(A, Ti[p-1 for p in perm]) end - else - return analyze_p(A, Ti[p-1 for p in perm]) end end end + end function cholesky!(F::Factor{Tv}, A::Sparse{Tv}; From ec38631e6b261eaea3a12be0764f18437e90c7a1 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" <744411+ViralBShah@users.noreply.github.com> Date: Sat, 2 Nov 2024 13:19:41 -0400 Subject: [PATCH 8/8] Update ci.yml with more architectures --- .github/workflows/ci.yml | 16 ++++------------ 1 file changed, 4 insertions(+), 12 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 4a6b9455..183f2ddc 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -28,12 +28,13 @@ jobs: - windows-latest arch: - x64 + - x86 include: - os: macOS-latest arch: aarch64 version: '1.10' - - os: ubuntu-latest - arch: x86 + - os: macOS-13 + arch: x64 version: '1.10' steps: - uses: actions/checkout@v4 @@ -41,16 +42,7 @@ jobs: with: version: ${{ matrix.version }} arch: ${{ matrix.arch }} - - uses: actions/cache@v4 - env: - cache-name: cache-artifacts - with: - path: ~/.julia/artifacts - key: ${{ runner.os }}-test-${{ env.cache-name }}-${{ hashFiles('**/Project.toml') }} - restore-keys: | - ${{ runner.os }}-test-${{ env.cache-name }}- - ${{ runner.os }}-test-${{ matrix.os }} - ${{ runner.os }}- + - uses: julia-actions/cache@v2 - run: julia --color=yes .ci/test_and_change_uuid.jl - uses: julia-actions/julia-buildpkg@v1 - uses: julia-actions/julia-runtest@v1