Skip to content

[WIP] Introduce AbstractSparseMatrixCSC into Type Hierarchy #24973

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

Closed
wants to merge 7 commits into from
Closed

[WIP] Introduce AbstractSparseMatrixCSC into Type Hierarchy #24973

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
ea7c32b
start converting SuiteSparse
JaredCrean2 Dec 6, 2017
4cfe084
workaround for type system invarience
JaredCrean2 Dec 7, 2017
af96299
misc fixes
JaredCrean2 Dec 7, 2017
8978119
add cholmod tests
JaredCrean2 Dec 7, 2017
04cec41
add spqr and umfpack tests
JaredCrean2 Dec 7, 2017
51fb20c
update nnz() signature
JaredCrean2 Dec 7, 2017
7cbdbae
fix whitespace
JaredCrean2 Dec 7, 2017
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
1 change: 1 addition & 0 deletions base/exports.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1253,6 +1253,7 @@ export
SparseArrays,
AbstractSparseArray,
AbstractSparseMatrix,
AbstractSparseMatrixCSC,
AbstractSparseVector,
SparseMatrixCSC,
SparseVector,
Expand Down
8 changes: 4 additions & 4 deletions base/sparse/linalg.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ import Base.LinAlg: checksquare

## sparse matrix multiplication

function (*)(A::SparseMatrixCSC{TvA,TiA}, B::SparseMatrixCSC{TvB,TiB}) where {TvA,TiA,TvB,TiB}
function (*)(A::AbstractSparseMatrixCSC{TvA,TiA}, B::AbstractSparseMatrixCSC{TvB,TiB}) where {TvA,TiA,TvB,TiB}
(*)(sppromote(A, B)...)
end
for f in (:A_mul_Bt, :A_mul_Bc,
Expand All @@ -17,7 +17,7 @@ for f in (:A_mul_Bt, :A_mul_Bc,
end
end

function sppromote(A::SparseMatrixCSC{TvA,TiA}, B::SparseMatrixCSC{TvB,TiB}) where {TvA,TiA,TvB,TiB}
function sppromote(A::AbstractSparseMatrixCSC{TvA,TiA}, B::AbstractSparseMatrixCSC{TvB,TiB}) where {TvA,TiA,TvB,TiB}
Tv = promote_type(TvA, TvB)
Ti = promote_type(TiA, TiB)
A = convert(SparseMatrixCSC{Tv,Ti}, A)
Expand Down Expand Up @@ -106,7 +106,7 @@ end
# Sparse matrix multiplication as described in [Gustavson, 1978]:
# http://dl.acm.org/citation.cfm?id=355796

(*)(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti}) where {Tv,Ti} = spmatmul(A,B)
(*)(A::AbstractSparseMatrixCSC{Tv,Ti}, B::AbstractSparseMatrixCSC{Tv,Ti}) where {Tv,Ti} = spmatmul(A,B)
for (f, opA, opB) in ((:A_mul_Bt, :identity, :transpose),
(:A_mul_Bc, :identity, :adjoint),
(:At_mul_B, :transpose, :identity),
Expand All @@ -120,7 +120,7 @@ for (f, opA, opB) in ((:A_mul_Bt, :identity, :transpose),
end
end

function spmatmul(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti};
function spmatmul(A::AbstractSparseMatrixCSC{Tv,Ti}, B::AbstractSparseMatrixCSC{Tv,Ti};
sortindices::Symbol = :sortcols) where {Tv,Ti}
mA, nA = size(A)
mB, nB = size(B)
Expand Down
3 changes: 2 additions & 1 deletion base/sparse/sparse.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,8 @@ import Base: @get!, acos, acosd, acot, acotd, acsch, asech, asin, asind, asinh,
rotl90, rotr90, round, scale!, setindex!, similar, size, transpose, tril,
triu, vec, permute!, map, map!

export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseVector,
export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseMatrixCSC,
AbstractSparseVector,
SparseMatrixCSC, SparseVector, blkdiag, droptol!, dropzeros!, dropzeros,
issparse, nonzeros, nzrange, rowvals, sparse, sparsevec, spdiagm, spones,
sprand, sprandn, spzeros, nnz, permute
Expand Down
13 changes: 11 additions & 2 deletions base/sparse/sparsematrix.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,13 +5,22 @@
# Assumes that row values in rowval for each column are sorted
# issorted(rowval[colptr[i]:(colptr[i+1]-1)]) == true

"""
AbstractSparseMatrixCSC{Tv,Ti} <: AbstractSparseMatrix{Tv,Ti}

Supertype of all sparse matrices that store the data using
[Compressed Sparse Column](@ref man-csc) format. Concrete types must have all the fields
of [`SparseMatrixCSC`](@ref) to be usable by external linear algebra libraries.
"""
abstract type AbstractSparseMatrixCSC{Tv,Ti} <: AbstractSparseMatrix{Tv,Ti} end

"""
SparseMatrixCSC{Tv,Ti<:Integer} <: AbstractSparseMatrix{Tv,Ti}

Matrix type for storing sparse matrices in the
[Compressed Sparse Column](@ref man-csc) format.
"""
struct SparseMatrixCSC{Tv,Ti<:Integer} <: AbstractSparseMatrix{Tv,Ti}
struct SparseMatrixCSC{Tv,Ti<:Integer} <: AbstractSparseMatrixCSC{Tv,Ti}
m::Int # Number of rows
n::Int # Number of columns
colptr::Vector{Ti} # Column i is in colptr[i]:(colptr[i+1]-1)
Expand Down Expand Up @@ -69,7 +78,7 @@ julia> nnz(A)
3
```
"""
nnz(S::SparseMatrixCSC) = Int(S.colptr[S.n + 1] - 1)
nnz(S::AbstractSparseMatrixCSC) = Int(S.colptr[S.n + 1] - 1)
count(pred, S::SparseMatrixCSC) = count(pred, view(S.nzval, 1:nnz(S))) + pred(zero(eltype(S)))*(prod(size(S)) - nnz(S))

"""
Expand Down
140 changes: 81 additions & 59 deletions stdlib/SuiteSparse/src/cholmod.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,9 @@ export
Factor,
Sparse

import ..SparseArrays: AbstractSparseMatrix, SparseMatrixCSC, indtype, sparse, spzeros, nnz
# why are these imports needed? most of these names are exported?
import ..SparseArrays: AbstractSparseMatrix, SparseMatrixCSC, AbstractSparseMatrixCSC,
indtype, sparse, spzeros, nnz

import ..increment, ..increment!, ..decrement, ..decrement!

Expand Down Expand Up @@ -249,6 +251,8 @@ struct C_SparseVoid <: SuiteSparseStruct
packed::Cint
end

# although this contains a SparseMatrixCSC via the C_SparseVoid struct, making
# it an AbstractSparseMatrixCSC changes method specificity in undesirable ways
mutable struct Sparse{Tv<:VTypes} <: AbstractSparseMatrix{Tv,SuiteSparse_long}
p::Ptr{C_Sparse{Tv}}
function Sparse{Tv}(p::Ptr{C_Sparse{Tv}}) where Tv<:VTypes
Expand Down Expand Up @@ -887,7 +891,7 @@ function Sparse(m::Integer, n::Integer,
o
end

function Sparse(A::SparseMatrixCSC{Tv,SuiteSparse_long}, stype::Integer) where Tv<:VTypes
function Sparse(A::AbstractSparseMatrixCSC{Tv,SuiteSparse_long}, stype::Integer) where Tv<:VTypes
## Check length of input. This should never fail but see #20024
if length(A.colptr) <= A.n
throw(ArgumentError("length of colptr must be at least size(A,2) + 1 = $(A.n + 1) but was $(length(A.colptr))"))
Expand Down Expand Up @@ -917,23 +921,25 @@ end
# convert SparseVectors into CHOLMOD Sparse types through a mx1 CSC matrix
convert(::Type{Sparse}, A::SparseVector{<:VTypes,SuiteSparse_long}) =
convert(Sparse, convert(SparseMatrixCSC, A))
function convert(::Type{Sparse}, A::SparseMatrixCSC{<:VTypes,<:ITypes})
function convert(::Type{Sparse}, A::AbstractSparseMatrixCSC{<:VTypes,<:ITypes})
o = Sparse(A, 0)
# check if array is symmetric and change stype if it is
if ishermitian(o)
change_stype!(o, -1)
end
o
end
convert(::Type{Sparse}, A::SparseMatrixCSC{Complex{Float32},<:ITypes}) =
convert(::Type{Sparse}, A::AbstractSparseMatrixCSC{Complex{Float32},<:ITypes}) =
convert(Sparse, convert(SparseMatrixCSC{Complex{Float64},SuiteSparse_long}, A))
convert(::Type{Sparse}, A::Symmetric{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}}) =
Sparse(A.data, A.uplo == 'L' ? -1 : 1)
convert(::Type{Sparse}, A::Symmetric{Float64, Tcsc}) where {
Tcsc<:AbstractSparseMatrixCSC{Float64,SuiteSparse_long}
} = Sparse(A.data, A.uplo == 'L' ? -1 : 1)
# The convert method for Hermitian is very similar to the general convert method, but we need to
# remove any non real elements in the diagonal because, in contrast to BLAS/LAPACK these are not
# ignored by CHOLMOD. If even tiny imaginary parts are present CHOLMOD will fail with a non-positive
# definite/zero pivot error.
function convert(::Type{Sparse}, AH::Hermitian{Tv,SparseMatrixCSC{Tv,SuiteSparse_long}}) where {Tv<:VTypes}
function convert(::Type{Sparse}, AH::Hermitian{Tv,Tcsc}) where {
Tv<:VTypes,Tcsc<:AbstractSparseMatrixCSC{Tv,SuiteSparse_long}}
A = AH.data

# Here we allocate a Symmetric/Hermitian CHOLMOD.Sparse matrix so we only need to copy
Expand All @@ -958,17 +964,18 @@ function convert(::Type{Sparse}, AH::Hermitian{Tv,SparseMatrixCSC{Tv,SuiteSparse
return o
end
function convert(::Type{Sparse},
A::Union{SparseMatrixCSC{BigFloat,Ti},
Symmetric{BigFloat,SparseMatrixCSC{BigFloat,Ti}},
Hermitian{Complex{BigFloat},SparseMatrixCSC{Complex{BigFloat},Ti}}},
args...) where Ti<:ITypes
A::Union{Tcsc,
Symmetric{BigFloat, Tcsc},
Hermitian{Complex{BigFloat}, TCcsc}},
args...) where {Ti<:ITypes, Tcsc<:AbstractSparseMatrixCSC{BigFloat,Ti},
TCcsc<:AbstractSparseMatrixCSC{Complex{BigFloat},Ti}}
throw(MethodError(convert, (Sparse, A)))
end
function convert(::Type{Sparse},
A::Union{SparseMatrixCSC{T,Ti},
Symmetric{T,SparseMatrixCSC{T,Ti}},
Hermitian{T,SparseMatrixCSC{T,Ti}}},
args...) where T where Ti<:ITypes
A::Union{Tcsc,
Symmetric{T,Tcsc},
Hermitian{T,Tcsc}},
args...) where {T, Ti<:ITypes,Tcsc<:AbstractSparseMatrixCSC{T,Ti}}
return Sparse(convert(AbstractMatrix{promote_type(Float64, T)}, A), args...)
end

Expand Down Expand Up @@ -1272,7 +1279,7 @@ function getindex(F::Factor, sym::Symbol)
FactorComponent(F, sym)
end

function getLd!(S::SparseMatrixCSC)
function getLd!(S::AbstractSparseMatrixCSC)
d = Vector{eltype(S)}(uninitialized, size(S, 1))
fill!(d, 0)
col = 1
Expand Down Expand Up @@ -1365,8 +1372,8 @@ end
cholfact!(F::Factor, A; shift = 0.0) -> CHOLMOD.Factor

Compute the Cholesky (``LL'``) factorization of `A`, reusing the symbolic
factorization `F`. `A` must be a [`SparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/
[`Hermitian`](@ref) view of a `SparseMatrixCSC`. Note that even if `A` doesn't
factorization `F`. `A` must be an [`AbstractSparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/
[`Hermitian`](@ref) view of a `AbstractSparseMatrixCSC`. Note that even if `A` doesn't
have the type tag, it must still be symmetric or Hermitian.

See also [`cholfact`](@ref).
Expand All @@ -1377,12 +1384,14 @@ See also [`cholfact`](@ref).
be converted to `SparseMatrixCSC{Float64}` or `SparseMatrixCSC{Complex128}`
as appropriate.
"""
cholfact!(F::Factor, A::Union{SparseMatrixCSC{T},
SparseMatrixCSC{Complex{T}},
Symmetric{T,SparseMatrixCSC{T,SuiteSparse_long}},
Hermitian{Complex{T},SparseMatrixCSC{Complex{T},SuiteSparse_long}},
Hermitian{T,SparseMatrixCSC{T,SuiteSparse_long}}};
shift = 0.0) where {T<:Real} =
cholfact!(F::Factor, A::Union{AbstractSparseMatrixCSC{T},
AbstractSparseMatrixCSC{Complex{T}},
Symmetric{T,Tcsc},
Hermitian{Complex{T},TCcsc},
Hermitian{T,Tcsc}};
shift = 0.0) where {T<:Real,
Tcsc<:AbstractSparseMatrixCSC{T,SuiteSparse_long},
TCcsc<:AbstractSparseMatrixCSC{Complex{T},SuiteSparse_long}} =
cholfact!(F, Sparse(A); shift = shift)

function cholfact(A::Sparse; shift::Real=0.0,
Expand All @@ -1404,8 +1413,8 @@ end
cholfact(A; shift = 0.0, perm = Int[]) -> CHOLMOD.Factor

Compute the Cholesky factorization of a sparse positive definite matrix `A`.
`A` must be a [`SparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `SparseMatrixCSC`. Note that even if `A` doesn't
`A` must be a [`AbstractSparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `AbstractSparseMatrixCSC`. Note that even if `A` doesn't
have the type tag, it must still be symmetric or Hermitian.
A fill-reducing permutation is used.
`F = cholfact(A)` is most frequently used to solve systems of equations with `F\\b`,
Expand Down Expand Up @@ -1433,11 +1442,13 @@ it should be a permutation of `1:size(A,1)` giving the ordering to use
Many other functions from CHOLMOD are wrapped but not exported from the
`Base.SparseArrays.CHOLMOD` module.
"""
cholfact(A::Union{SparseMatrixCSC{T}, SparseMatrixCSC{Complex{T}},
Symmetric{T,SparseMatrixCSC{T,SuiteSparse_long}},
Hermitian{Complex{T},SparseMatrixCSC{Complex{T},SuiteSparse_long}},
Hermitian{T,SparseMatrixCSC{T,SuiteSparse_long}}};
kws...) where {T<:Real} = cholfact(Sparse(A); kws...)
cholfact(A::Union{AbstractSparseMatrixCSC{T}, AbstractSparseMatrixCSC{Complex{T}},
Symmetric{T, Tcsc},
Hermitian{Complex{T},TCcsc},
Hermitian{T,Tcsc}};
kws...) where {T<:Real, Tcsc<:AbstractSparseMatrixCSC{T, SuiteSparse_long},
TCcsc<:AbstractSparseMatrixCSC{Complex{T}, SuiteSparse_long}
} = cholfact(Sparse(A); kws...)


function ldltfact!(F::Factor{Tv}, A::Sparse{Tv}; shift::Real=0.0) where Tv
Expand All @@ -1457,8 +1468,8 @@ end
ldltfact!(F::Factor, A; shift = 0.0) -> CHOLMOD.Factor

Compute the ``LDL'`` factorization of `A`, reusing the symbolic factorization `F`.
`A` must be a [`SparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `SparseMatrixCSC`. Note that even if `A` doesn't
`A` must be a [`AbstractSparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `AbstractSparseMatrixCSC`. Note that even if `A` doesn't
have the type tag, it must still be symmetric or Hermitian.

See also [`ldltfact`](@ref).
Expand All @@ -1469,12 +1480,14 @@ See also [`ldltfact`](@ref).
be converted to `SparseMatrixCSC{Float64}` or `SparseMatrixCSC{Complex128}`
as appropriate.
"""
ldltfact!(F::Factor, A::Union{SparseMatrixCSC{T},
SparseMatrixCSC{Complex{T}},
Symmetric{T,SparseMatrixCSC{T,SuiteSparse_long}},
Hermitian{Complex{T},SparseMatrixCSC{Complex{T},SuiteSparse_long}},
Hermitian{T,SparseMatrixCSC{T,SuiteSparse_long}}};
shift = 0.0) where {T<:Real} =
ldltfact!(F::Factor, A::Union{AbstractSparseMatrixCSC{T},
AbstractSparseMatrixCSC{Complex{T}},
Symmetric{T,Tcsc},
Hermitian{Complex{T},TCcsc},
Hermitian{T,Tcsc}};
shift = 0.0) where {T<:Real,
Tcsc<:AbstractSparseMatrixCSC{T,SuiteSparse_long},
TCcsc<:AbstractSparseMatrixCSC{Complex{T},SuiteSparse_long}} =
ldltfact!(F, Sparse(A), shift = shift)

function ldltfact(A::Sparse; shift::Real=0.0,
Expand All @@ -1501,8 +1514,8 @@ end
ldltfact(A; shift = 0.0, perm=Int[]) -> CHOLMOD.Factor

Compute the ``LDL'`` factorization of a sparse matrix `A`.
`A` must be a [`SparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `SparseMatrixCSC`. Note that even if `A` doesn't
`A` must be a [`AbstractSparseMatrixCSC`](@ref) or a [`Symmetric`](@ref)/[`Hermitian`](@ref)
view of a `AbstractSparseMatrixCSC`. Note that even if `A` doesn't
have the type tag, it must still be symmetric or Hermitian.
A fill-reducing permutation is used. `F = ldltfact(A)` is most frequently
used to solve systems of equations `A*x = b` with `F\\b`. The returned
Expand Down Expand Up @@ -1531,11 +1544,14 @@ it should be a permutation of `1:size(A,1)` giving the ordering to use
Many other functions from CHOLMOD are wrapped but not exported from the
`Base.SparseArrays.CHOLMOD` module.
"""
ldltfact(A::Union{SparseMatrixCSC{T},SparseMatrixCSC{Complex{T}},
Symmetric{T,SparseMatrixCSC{T,SuiteSparse_long}},
Hermitian{Complex{T},SparseMatrixCSC{Complex{T},SuiteSparse_long}},
Hermitian{T,SparseMatrixCSC{T,SuiteSparse_long}}};
kws...) where {T<:Real} = ldltfact(Sparse(A); kws...)
ldltfact(A::Union{AbstractSparseMatrixCSC{T},AbstractSparseMatrixCSC{Complex{T}},
Symmetric{T,Tcsc},
Hermitian{Complex{T},TCcsc},
Hermitian{T,Tcsc}};
kws...) where {T<:Real,
Tcsc<:AbstractSparseMatrixCSC{T,SuiteSparse_long},
TCcsc<:AbstractSparseMatrixCSC{Complex{T},SuiteSparse_long}} =
ldltfact(Sparse(A); kws...)

## Rank updates

Expand Down Expand Up @@ -1660,7 +1676,7 @@ for (T, f) in ((:Dense, :solve), (:Sparse, :spsolve))
end
end

SparseVecOrMat{Tv,Ti} = Union{SparseVector{Tv,Ti}, SparseMatrixCSC{Tv,Ti}}
SparseVecOrMat{Tv,Ti} = Union{SparseVector{Tv,Ti}, AbstractSparseMatrixCSC{Tv,Ti}}

function (\)(L::FactorComponent, b::Vector)
reshape(convert(Matrix, L\Dense(b)), length(b))
Expand Down Expand Up @@ -1692,9 +1708,12 @@ Ac_ldiv_B(L::Factor, B::Sparse) = spsolve(CHOLMOD_A, L, B)
Ac_ldiv_B(L::Factor, B::SparseVecOrMat) = Ac_ldiv_B(L, Sparse(B))

for f in (:\, :Ac_ldiv_B)
@eval function ($f)(A::Union{Symmetric{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}},
Hermitian{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}},
Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},SuiteSparse_long}}}, B::StridedVecOrMat)
@eval function ($f)(A::Union{Symmetric{Float64, Tcsc},
Hermitian{Float64,Tcsc},
Hermitian{Complex{Float64}, TCcsc}},
B::StridedVecOrMat) where {
Tcsc<:AbstractSparseMatrixCSC{Float64,SuiteSparse_long},
TCcsc<:AbstractSparseMatrixCSC{Complex{Float64},SuiteSparse_long}}
F = cholfact(A)
if issuccess(F)
return ($f)(F, B)
Expand Down Expand Up @@ -1801,18 +1820,21 @@ function ishermitian(A::Sparse{Complex{Float64}})
end
end

(*)(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}},
B::SparseVecOrMat{Float64,Ti}) where {Ti} = sparse(Sparse(A)*Sparse(B))
(*)(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}},
B::SparseVecOrMat{Complex{Float64},Ti}) where {Ti} = sparse(Sparse(A)*Sparse(B))
(*)(A::Hermitian{Float64,SparseMatrixCSC{Float64,Ti}},
B::SparseVecOrMat{Float64,Ti}) where {Ti} = sparse(Sparse(A)*Sparse(B))
#TODO: make these return in same type as A (which might not be SparseMatrixCSC).
(*)(A::Symmetric{Float64, Tcsc},
B::SparseVecOrMat{Float64,Ti}) where {
Ti, Tcsc<:AbstractSparseMatrixCSC{Float64,Ti}} = sparse(Sparse(A)*Sparse(B))
(*)(A::Hermitian{Complex{Float64},TCcsc},
B::SparseVecOrMat{Complex{Float64},Ti}) where {Ti,
TCcsc<:AbstractSparseMatrixCSC{Complex{Float64},Ti}} = sparse(Sparse(A)*Sparse(B))
(*)(A::Hermitian{Float64,Tcsc},
B::SparseVecOrMat{Float64,Ti}) where {Ti, Tcsc<:AbstractSparseMatrixCSC{Float64,Ti}} = sparse(Sparse(A)*Sparse(B))

(*)(A::SparseVecOrMat{Float64,Ti},
B::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}) where {Ti} = sparse(Sparse(A)*Sparse(B))
B::Symmetric{Float64,Tcsc}) where {Ti, Tcsc<:AbstractSparseMatrixCSC{Float64,Ti}} = sparse(Sparse(A)*Sparse(B))
(*)(A::SparseVecOrMat{Complex{Float64},Ti},
B::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}) where {Ti} = sparse(Sparse(A)*Sparse(B))
B::Hermitian{Complex{Float64},TCcsc}) where {Ti, TCcsc<:AbstractSparseMatrixCSC{Complex{Float64},Ti}} = sparse(Sparse(A)*Sparse(B))
(*)(A::SparseVecOrMat{Float64,Ti},
B::Hermitian{Float64,SparseMatrixCSC{Float64,Ti}}) where {Ti} = sparse(Sparse(A)*Sparse(B))
B::Hermitian{Float64,Tcsc}) where {Ti, Tcsc<:AbstractSparseMatrixCSC{Float64,Ti}} = sparse(Sparse(A)*Sparse(B))

end #module
10 changes: 5 additions & 5 deletions stdlib/SuiteSparse/src/spqr.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ const ORDERING_BESTAMD = Int32(9) # try COLAMD and AMD; pick best#
# tried. If there is a high fill-in with AMD then try METIS(A'A) and take
# the best of AMD and METIS. METIS is not tried if it isn't installed.

using ..SparseArrays: SparseMatrixCSC
using ..SparseArrays: SparseMatrixCSC, AbstractSparseMatrixCSC
using ..SuiteSparse.CHOLMOD
using ..SuiteSparse.CHOLMOD: change_stype!, free!

Expand Down Expand Up @@ -132,10 +132,10 @@ end
Base.size(Q::QRSparseQ) = (size(Q.factors, 1), size(Q.factors, 1))

# From SPQR manual p. 6
_default_tol(A::SparseMatrixCSC) =
_default_tol(A::AbstractSparseMatrixCSC) =
20*sum(size(A))*eps(real(eltype(A)))*maximum(norm(view(A, :, i))^2 for i in 1:size(A, 2))

function Base.LinAlg.qrfact(A::SparseMatrixCSC{Tv}; tol = _default_tol(A)) where {Tv <: CHOLMOD.VTypes}
function Base.LinAlg.qrfact(A::AbstractSparseMatrixCSC{Tv}; tol = _default_tol(A)) where {Tv <: CHOLMOD.VTypes}
R = Ref{Ptr{CHOLMOD.C_Sparse{Tv}}}()
E = Ref{Ptr{CHOLMOD.SuiteSparse_long}}()
H = Ref{Ptr{CHOLMOD.C_Sparse{Tv}}}()
Expand Down Expand Up @@ -193,9 +193,9 @@ Column permutation:
2
```
"""
Base.LinAlg.qrfact(A::SparseMatrixCSC; tol = _default_tol(A)) = qrfact(A, Val{true}, tol = tol)
Base.LinAlg.qrfact(A::AbstractSparseMatrixCSC; tol = _default_tol(A)) = qrfact(A, Val{true}, tol = tol)

Base.LinAlg.qr(A::SparseMatrixCSC; tol = _default_tol(A)) = qr(A, Val{true}, tol = tol)
Base.LinAlg.qr(A::AbstractSparseMatrixCSC; tol = _default_tol(A)) = qr(A, Val{true}, tol = tol)

function Base.A_mul_B!(Q::QRSparseQ, A::StridedVecOrMat)
if size(A, 1) != size(Q, 1)
Expand Down
Loading