Skip to content

Commit

Permalink
Merge pull request #439 from JuliaSparse/backports-release-1.10
Browse files Browse the repository at this point in the history
Backports for Julia v1.10
  • Loading branch information
dkarrasch authored Oct 11, 2023
2 parents 026e6a6 + a272109 commit d8aae6d
Show file tree
Hide file tree
Showing 10 changed files with 114 additions and 51 deletions.
3 changes: 3 additions & 0 deletions docs/Project.toml
Original file line number Diff line number Diff line change
@@ -1,2 +1,5 @@
[deps]
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"

[compat]
Documenter = "1"
3 changes: 1 addition & 2 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,7 @@ makedocs(
"SparseArrays" => "index.md",
"Sparse Linear Algebra" => "solvers.md",
];
# strict = true,
strict = Symbol[:doctest],
warnonly = [:missing_docs, :cross_references],
)

deploydocs(repo = "github.com/JuliaSparse/SparseArrays.jl.git")
13 changes: 8 additions & 5 deletions docs/src/solvers.md
Original file line number Diff line number Diff line change
Expand Up @@ -14,11 +14,14 @@ Sparse matrix solvers call functions from [SuiteSparse](http://suitesparse.com).

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 [`Linear Algebra`](https://docs.julialang.org/en/v1/stdlib/LinearAlgebra/) section of the manual:
1. [`cholesky`](@ref)
2. [`ldlt`](@ref)
3. [`lu`](@ref)
4. [`qr`](@ref)
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
SparseArrays.CHOLMOD.cholesky
Expand Down
2 changes: 1 addition & 1 deletion src/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -560,7 +560,7 @@ function LinearAlgebra.generic_trimatmul!(C::StridedVecOrMat, uploc, isunitc, tf
end
return C
end
function LinearAlgebra.generic_trimatmul!(C::StridedVecOrMat, uploc, isunitc, _, xA::AdjOrTrans{<:Any,<:SparseMatrixCSCUnion}, B::AbstractVecOrMat)
function LinearAlgebra.generic_trimatmul!(C::StridedVecOrMat, uploc, isunitc, ::Function, xA::AdjOrTrans{<:Any,<:SparseMatrixCSCUnion}, B::AbstractVecOrMat)
A = parent(xA)
nrowC = size(C, 1)
ncol = checksquare(A)
Expand Down
7 changes: 4 additions & 3 deletions src/solvers/cholmod.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,8 @@ using LinearAlgebra
using LinearAlgebra: RealHermSymComplexHerm, AdjOrTrans
import LinearAlgebra: (\), AdjointFactorization,
cholesky, cholesky!, det, diag, ishermitian, isposdef,
issuccess, issymmetric, ldlt, ldlt!, logdet, lowrankdowndate!
issuccess, issymmetric, ldlt, ldlt!, logdet,
lowrankdowndate, lowrankdowndate!, lowrankupdate, lowrankupdate!

using SparseArrays
using SparseArrays: getcolptr, AbstractSparseVecOrMat
Expand Down Expand Up @@ -1549,7 +1550,7 @@ factor will be `L*L' == P*A*P' + C'*C`
`update`: `Cint(1)` for `A + CC'`, `Cint(0)` for `A - CC'`
"""
lowrankdowndate!
lowrankupdowndate!

#Helper functions for rank updates
lowrank_reorder(V::AbstractArray,p) = Sparse(sparse(V[p,:]))
Expand Down Expand Up @@ -1598,7 +1599,7 @@ lowrankupdate(F::Factor{Tv}, V::AbstractArray{Tv}) where {Tv<:VTypes} =
lowrankupdate!(copy(F), V)

"""
lowrankupdate(F::CHOLMOD.Factor, C::AbstractArray) -> FF::CHOLMOD.Factor
lowrankdowndate(F::CHOLMOD.Factor, C::AbstractArray) -> FF::CHOLMOD.Factor
Get an `LDLt` Factorization of `A + C*C'` given an `LDLt` or `LLt` factorization `F` of `A`.
Expand Down
36 changes: 21 additions & 15 deletions src/sparsematrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -344,6 +344,26 @@ function Base.print_array(io::IO, S::AbstractSparseMatrixCSCInclAdjointAndTransp
end
end

"""
ColumnIndices(S::AbstractSparseMatrixCSC)
Return the column indices of the stored values in `S`.
This is an internal type that is used in displaying sparse matrices,
and is not a part of the public interface.
"""
struct ColumnIndices{Ti,S<:AbstractSparseMatrixCSC{<:Any,Ti}} <: AbstractVector{Ti}
arr :: S
end

size(C::ColumnIndices) = (nnz(C.arr),)
# returns the column index of the n-th non-zero value from the column pointer
@inline function getindex(C::ColumnIndices, i::Int)
@boundscheck checkbounds(C, i)
colptr = getcolptr(C.arr)
ind = searchsortedlast(colptr, i)
eltype(C)(ind)
end

# always show matrices as `sparse(I, J, K)`
function Base.show(io::IO, _S::AbstractSparseMatrixCSCInclAdjointAndTranspose)
_checkbuffers(_S)
Expand All @@ -358,21 +378,7 @@ function Base.show(io::IO, _S::AbstractSparseMatrixCSCInclAdjointAndTranspose)
print(io, "transpose(")
end
print(io, "sparse(", I, ", ")
if length(I) == 0
print(io, eltype(getcolptr(S)), "[]")
else
print(io, "[")
il = nnz(S) - 1
for col in 1:size(S, 2),
k in getcolptr(S)[col] : (getcolptr(S)[col+1]-1)
print(io, col)
if il > 0
print(io, ", ")
il -= 1
end
end
print(io, "]")
end
show(io, ColumnIndices(S))
print(io, ", ", K, ", ", m, ", ", n, ")")
if _S isa Adjoint || _S isa Transpose
print(io, ")")
Expand Down
57 changes: 36 additions & 21 deletions src/sparsevector.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1108,7 +1108,9 @@ function hcat(Xin::AbstractSparseVector...)
X = map(_unsafe_unfix, Xin)
Tv = promote_type(map(eltype, X)...)
Ti = promote_type(map(indtype, X)...)
r = _absspvec_hcat(map(x -> convert(SparseVector{Tv,Ti}, x), X)...)
r = (function (::Type{SV}) where SV
_absspvec_hcat(map(x -> convert(SV, x), X)...)
end)(SparseVector{Tv,Ti})
return @if_move_fixed Xin... r
end
function _absspvec_hcat(X::AbstractSparseVector{Tv,Ti}...) where {Tv,Ti}
Expand Down Expand Up @@ -1144,7 +1146,9 @@ function vcat(Xin::AbstractSparseVector...)
X = map(_unsafe_unfix, Xin)
Tv = promote_type(map(eltype, X)...)
Ti = promote_type(map(indtype, X)...)
r = _absspvec_vcat(map(x -> convert(SparseVector{Tv,Ti}, x), X)...)
r = (function (::Type{SV}) where SV
_absspvec_vcat(map(x -> convert(SV, x), X)...)
end)(SparseVector{Tv,Ti})
return @if_move_fixed Xin... r
end
function _absspvec_vcat(X::AbstractSparseVector{Tv,Ti}...) where {Tv,Ti}
Expand Down Expand Up @@ -1193,14 +1197,16 @@ _makesparse(x::AbstractMatrix) = convert(SparseMatrixCSC, issparse(x) ? x : spar
anysparse() = false
anysparse(X) = X isa AbstractArray && issparse(X)
anysparse(X, Xs...) = anysparse(X) || anysparse(Xs...)
anysparse(X::T, Xs::T...) where {T} = anysparse(X)

function hcat(X::Union{Vector, AbstractSparseVector}...)
const _SparseVecConcatGroup = Union{Vector, AbstractSparseVector}
function hcat(X::_SparseVecConcatGroup...)
if anysparse(X...)
X = map(sparse, X)
end
return cat(X...; dims=Val(2))
end
function vcat(X::Union{Vector, AbstractSparseVector}...)
function vcat(X::_SparseVecConcatGroup...)
if anysparse(X...)
X = map(sparse, X)
end
Expand All @@ -1213,30 +1219,30 @@ end
const _SparseConcatGroup = Union{AbstractVecOrMat{<:Number},Number}

# `@constprop :aggressive` allows `dims` to be propagated as constant improving return type inference
Base.@constprop :aggressive function Base._cat(dims, X::_SparseConcatGroup...)
T = promote_eltype(X...)
if anysparse(X...)
X = (_sparse(first(X)), map(_makesparse, Base.tail(X))...)
Base.@constprop :aggressive function Base._cat(dims, X1::_SparseConcatGroup, X::_SparseConcatGroup...)
T = promote_eltype(X1, X...)
if anysparse(X1) || anysparse(X...)
X1, X = _sparse(X1), map(_makesparse, X)
end
return Base._cat_t(dims, T, X...)
return Base._cat_t(dims, T, X1, X...)
end
function hcat(X::_SparseConcatGroup...)
if anysparse(X...)
X = (_sparse(first(X)), map(_makesparse, Base.tail(X))...)
function hcat(X1::_SparseConcatGroup, X::_SparseConcatGroup...)
if anysparse(X1) || anysparse(X...)
X1, X = _sparse(X1), map(_makesparse, X)
end
return cat(X..., dims=Val(2))
return Base.typed_hcat(Base.promote_eltype(X1, X...), X1, X...)
end
function vcat(X::_SparseConcatGroup...)
if anysparse(X...)
X = (_sparse(first(X)), map(_makesparse, Base.tail(X))...)
function vcat(X1::_SparseConcatGroup, X::_SparseConcatGroup...)
if anysparse(X1) || anysparse(X...)
X1, X = _sparse(X1), map(_makesparse, X)
end
return cat(X..., dims=Val(1))
return Base.typed_vcat(Base.promote_eltype(X1, X...), X1, X...)
end
function hvcat(rows::Tuple{Vararg{Int}}, X::_SparseConcatGroup...)
if anysparse(X...)
vcat(_hvcat_rows(rows, X...)...)
function hvcat(rows::Tuple{Vararg{Int}}, X1::_SparseConcatGroup, X::_SparseConcatGroup...)
if anysparse(X1) || anysparse(X...)
vcat(_hvcat_rows(rows, X1, X...)...)
else
Base.typed_hvcat(Base.promote_eltypeof(X...), rows, X...)
Base.typed_hvcat(Base.promote_eltypeof(X1, X...), rows, X1, X...)
end
end
function _hvcat_rows((row1, rows...)::Tuple{Vararg{Int}}, X::_SparseConcatGroup...)
Expand All @@ -1254,6 +1260,15 @@ function _hvcat_rows((row1, rows...)::Tuple{Vararg{Int}}, X::_SparseConcatGroup.
end
_hvcat_rows(::Tuple{}, X::_SparseConcatGroup...) = ()

# disambiguation for type-piracy problems created above
hcat(n1::Number, ns::Vararg{Number}) = invoke(hcat, Tuple{Vararg{Number}}, n1, ns...)
vcat(n1::Number, ns::Vararg{Number}) = invoke(vcat, Tuple{Vararg{Number}}, n1, ns...)
hcat(n1::Type{N}, ns::Vararg{N}) where {N<:Number} = invoke(hcat, Tuple{Vararg{Number}}, n1, ns...)
vcat(n1::Type{N}, ns::Vararg{N}) where {N<:Number} = invoke(vcat, Tuple{Vararg{Number}}, n1, ns...)
hvcat(rows::Tuple{Vararg{Int}}, n1::Number, ns::Vararg{Number}) = invoke(hvcat, Tuple{typeof(rows), Vararg{Number}}, rows, n1, ns...)
hvcat(rows::Tuple{Vararg{Int}}, n1::N, ns::Vararg{N}) where {N<:Number} = invoke(hvcat, Tuple{typeof(rows), Vararg{N}}, rows, n1, ns...)


# make sure UniformScaling objects are converted to sparse matrices for concatenation
promote_to_array_type(A::Tuple{Vararg{Union{_SparseConcatGroup,UniformScaling}}}) = anysparse(A...) ? SparseMatrixCSC : Matrix
promote_to_arrays_(n::Int, ::Type{SparseMatrixCSC}, J::UniformScaling) = sparse(J, n, n)
Expand Down
23 changes: 19 additions & 4 deletions test/ambiguous.jl
Original file line number Diff line number Diff line change
Expand Up @@ -25,13 +25,13 @@ using Test, LinearAlgebra, SparseArrays, Aqua

@testset "code quality" begin
@testset "Method ambiguity" begin
# Aqua.test_ambiguities([SparseArrays, Base, Core])
Aqua.test_ambiguities([SparseArrays, Base, Core])
end
@testset "Unbound type parameters" begin
@test_broken Aqua.detect_unbound_args_recursively(SparseArrays) == []
end
@testset "Undefined exports" begin
Aqua.test_undefined_exports(SparseArrays) == []
Aqua.test_undefined_exports(SparseArrays)
end
@testset "Compare Project.toml and test/Project.toml" begin
Aqua.test_project_extras(SparseArrays)
Expand All @@ -50,8 +50,23 @@ using Test, LinearAlgebra, SparseArrays, Aqua
end
end

@testset "detect_ambiguities" begin
@test_nowarn detect_ambiguities(SparseArrays; recursive=true, ambiguous_bottom=false)
let ambig = detect_ambiguities(SparseArrays; recursive=true)
@test isempty(ambig)
ambig = Set{Any}(((m1.sig, m2.sig) for (m1, m2) in ambig))
expect = []
good = true
while !isempty(ambig)
sigs = pop!(ambig)
i = findfirst(==(sigs), expect)
if i === nothing
println(stderr, "push!(expect, (", sigs[1], ", ", sigs[2], "))")
good = false
continue
end
deleteat!(expect, i)
end
@test isempty(expect)
@test good
end

## This was the older version that was disabled
Expand Down
12 changes: 12 additions & 0 deletions test/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -162,6 +162,18 @@ begin
vMAW = tr(wr(view(MA, :, 1:n)))
@test vAW * B vMAW * B
end
a = sprand(rng, ComplexF64, n, n, 0.01)
ma = Matrix(a)
@testset "triangular multiply with conjugate matrices" for tr in (x -> adjoint(transpose(x)), x -> transpose(adjoint(x))),
wr in (UpperTriangular, LowerTriangular, UnitUpperTriangular, UnitLowerTriangular)
AW = tr(wr(a))
MAW = tr(wr(ma))
@test AW * B MAW * B
# and for SparseMatrixCSCView - a view of all rows and unit range of cols
vAW = tr(wr(view(a, :, 1:n)))
vMAW = tr(wr(view(ma, :, 1:n)))
@test vAW * B vMAW * B
end
A = A - Diagonal(diag(A)) + 2I # avoid rounding errors by division
MA = Matrix(A)
@testset "triangular solver for $tr($wr)" for tr in (identity, adjoint, transpose),
Expand Down
9 changes: 9 additions & 0 deletions test/sparsematrix_constructors_indexing.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1574,6 +1574,15 @@ end
_show_with_braille_patterns(ioc, _filled_sparse(8, 16))
@test String(take!(io)) == "⎡⣿⣿⎤\n" *
"⎣⣿⣿⎦"

# respect IOContext while displaying J
I, J, V = shuffle(1:50), shuffle(1:50), [1:50;]
S = sparse(I, J, V)
I, J, V = I[sortperm(J)], sort(J), V[sortperm(J)]
@test repr(S) == "sparse($I, $J, $V, $(size(S,1)), $(size(S,2)))"
limctxt(x) = repr(x, context=:limit=>true)
expstr = "sparse($(limctxt(I)), $(limctxt(J)), $(limctxt(V)), $(size(S,1)), $(size(S,2)))"
@test limctxt(S) == expstr
end

@testset "issparse for specialized matrix types" begin
Expand Down

0 comments on commit d8aae6d

Please sign in to comment.