Skip to content

Commit 760989e

Browse files
dkarraschdkarrasch
authored
Make sure transpose(UmfpackLU) is Transpose[Factorization] (#347)
1 parent 1c46709 commit 760989e

File tree

4 files changed

+29
-24
lines changed

4 files changed

+29
-24
lines changed

‎src/SparseArrays.jl

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -57,10 +57,10 @@ decrement(A::AbstractArray) = let y = Array(A)
5757
y .= y .- oneunit(eltype(A))
5858
end
5959

60-
AdjType = isdefined(LinearAlgebra, :AdjointFactorization) ?
60+
AdjointFact = isdefined(LinearAlgebra, :AdjointFactorization) ?
6161
LinearAlgebra.AdjointFactorization :
6262
Adjoint
63-
TransType = isdefined(LinearAlgebra, :TransposeFactorization) ?
63+
TransposeFact = isdefined(LinearAlgebra, :TransposeFactorization) ?
6464
LinearAlgebra.TransposeFactorization :
6565
Transpose
6666

‎src/solvers/cholmod.jl

Lines changed: 10 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -31,7 +31,7 @@ export
3131
import SparseArrays: AbstractSparseMatrix, SparseMatrixCSC, indtype, sparse, spzeros, nnz,
3232
sparsevec
3333

34-
import ..increment, ..increment!, ..AdjType, ..TransType
34+
import ..increment, ..increment!, ..AdjointFact, ..TransposeFact
3535

3636
using ..LibSuiteSparse
3737
import ..LibSuiteSparse: SuiteSparse_long, TRUE, FALSE
@@ -1577,21 +1577,21 @@ end
15771577
(\)(L::Factor, B::SparseVector) = sparsevec(spsolve(CHOLMOD_A, L, Sparse(B)))
15781578

15791579
# the eltype restriction is necessary for disambiguation with the B::StridedMatrix below
1580-
\(adjL::AdjType{<:VTypes,<:Factor}, B::Dense) = (L = adjL.parent; solve(CHOLMOD_A, L, B))
1581-
\(adjL::AdjType{<:Any,<:Factor}, B::Sparse) = (L = adjL.parent; spsolve(CHOLMOD_A, L, B))
1582-
\(adjL::AdjType{<:Any,<:Factor}, B::SparseVecOrMat) = (L = adjL.parent; \(adjoint(L), Sparse(B)))
1580+
\(adjL::AdjointFact{<:VTypes,<:Factor}, B::Dense) = (L = adjL.parent; solve(CHOLMOD_A, L, B))
1581+
\(adjL::AdjointFact{<:Any,<:Factor}, B::Sparse) = (L = adjL.parent; spsolve(CHOLMOD_A, L, B))
1582+
\(adjL::AdjointFact{<:Any,<:Factor}, B::SparseVecOrMat) = (L = adjL.parent; \(adjoint(L), Sparse(B)))
15831583

15841584
# Explicit typevars are necessary to avoid ambiguities with defs in LinearAlgebra/factorizations.jl
15851585
# Likewise the two following explicit Vector and Matrix defs (rather than a single VecOrMat)
1586-
(\)(adjL::AdjType{T,<:Factor}, B::Vector{Complex{T}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1587-
(\)(adjL::AdjType{T,<:Factor}, B::Matrix{Complex{T}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1588-
(\)(adjL::AdjType{T,<:Factor}, B::Adjoint{<:Any,Matrix{Complex{T}}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1589-
(\)(adjL::AdjType{T,<:Factor}, B::Transpose{<:Any,Matrix{Complex{T}}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1590-
function \(adjL::AdjType{<:VTypes,<:Factor}, b::StridedVector)
1586+
(\)(adjL::AdjointFact{T,<:Factor}, B::Vector{Complex{T}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1587+
(\)(adjL::AdjointFact{T,<:Factor}, B::Matrix{Complex{T}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1588+
(\)(adjL::AdjointFact{T,<:Factor}, B::Adjoint{<:Any,Matrix{Complex{T}}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1589+
(\)(adjL::AdjointFact{T,<:Factor}, B::Transpose{<:Any,Matrix{Complex{T}}}) where {T<:Float64} = complex.(adjL\real(B), adjL\imag(B))
1590+
function \(adjL::AdjointFact{<:VTypes,<:Factor}, b::StridedVector)
15911591
L = adjL.parent
15921592
return Vector(solve(CHOLMOD_A, L, Dense(b)))
15931593
end
1594-
function \(adjL::AdjType{<:VTypes,<:Factor}, B::StridedMatrix)
1594+
function \(adjL::AdjointFact{<:VTypes,<:Factor}, B::StridedMatrix)
15951595
L = adjL.parent
15961596
return Matrix(solve(CHOLMOD_A, L, Dense(B)))
15971597
end

‎src/solvers/umfpack.jl

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@ import SparseArrays: nnz
1515
import Serialization: AbstractSerializer, deserialize, serialize
1616
using Serialization
1717

18-
import ..increment, ..increment!, ..decrement, ..decrement!, ..AdjType, ..TransType
18+
import ..increment, ..increment!, ..decrement, ..decrement!, ..AdjointFact, ..TransposeFact
1919

2020
using ..LibSuiteSparse
2121
import ..LibSuiteSparse:
@@ -255,7 +255,7 @@ workspace_W_size(F::UmfpackLU) = workspace_W_size(F, has_refinement(F))
255255
workspace_W_size(S::Union{UmfpackLU{<:AbstractFloat}, AbstractSparseMatrixCSC{<:AbstractFloat}}, refinement::Bool) = refinement ? 5 * size(S, 2) : size(S, 2)
256256
workspace_W_size(S::Union{UmfpackLU{<:Complex}, AbstractSparseMatrixCSC{<:Complex}}, refinement::Bool) = refinement ? 10 * size(S, 2) : 4 * size(S, 2)
257257

258-
const ATLU = Union{TransType{<:Any, <:UmfpackLU}, AdjType{<:Any, <:UmfpackLU}}
258+
const ATLU = Union{TransposeFact{<:Any, <:UmfpackLU}, AdjointFact{<:Any, <:UmfpackLU}}
259259
has_refinement(F::ATLU) = has_refinement(F.parent)
260260
has_refinement(F::UmfpackLU) = has_refinement(F.control)
261261
has_refinement(control::AbstractVector) = control[JL_UMFPACK_IRSTEP] > 0
@@ -295,8 +295,8 @@ Base.copy(F::T, ws=UmfpackWS(F)) where {T <: ATLU} =
295295

296296
if !isdefined(LinearAlgebra, :AdjointFactorization)
297297
Base.adjoint(F::UmfpackLU) = Adjoint(F)
298-
Base.transpose(F::UmfpackLU) = Transpose(F)
299298
end
299+
Base.transpose(F::UmfpackLU) = TransposeFact(F)
300300

301301
function Base.lock(f::Function, F::UmfpackLU)
302302
lock(F)
@@ -938,28 +938,28 @@ import LinearAlgebra.ldiv!
938938

939939
ldiv!(lu::UmfpackLU{T}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
940940
ldiv!(B, lu, copy(B))
941-
ldiv!(translu::TransType{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
941+
ldiv!(translu::TransposeFact{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
942942
ldiv!(B, translu, copy(B))
943-
ldiv!(adjlu::AdjType{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
943+
ldiv!(adjlu::AdjointFact{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
944944
ldiv!(B, adjlu, copy(B))
945945
ldiv!(lu::UmfpackLU{Float64}, B::StridedVecOrMat{<:Complex}) =
946946
ldiv!(B, lu, copy(B))
947-
ldiv!(translu::TransType{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{<:Complex}) =
947+
ldiv!(translu::TransposeFact{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{<:Complex}) =
948948
ldiv!(B, translu, copy(B))
949-
ldiv!(adjlu::AdjType{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{<:Complex}) =
949+
ldiv!(adjlu::AdjointFact{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{<:Complex}) =
950950
ldiv!(B, adjlu, copy(B))
951951

952952
ldiv!(X::StridedVecOrMat{T}, lu::UmfpackLU{T}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
953953
_Aq_ldiv_B!(X, lu, B, UMFPACK_A)
954-
ldiv!(X::StridedVecOrMat{T}, translu::TransType{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
954+
ldiv!(X::StridedVecOrMat{T}, translu::TransposeFact{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
955955
(lu = translu.parent; _Aq_ldiv_B!(X, lu, B, UMFPACK_Aat))
956-
ldiv!(X::StridedVecOrMat{T}, adjlu::AdjType{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
956+
ldiv!(X::StridedVecOrMat{T}, adjlu::AdjointFact{T,<:UmfpackLU{T}}, B::StridedVecOrMat{T}) where {T<:UMFVTypes} =
957957
(lu = adjlu.parent; _Aq_ldiv_B!(X, lu, B, UMFPACK_At))
958958
ldiv!(X::StridedVecOrMat{Tb}, lu::UmfpackLU{Float64}, B::StridedVecOrMat{Tb}) where {Tb<:Complex} =
959959
_Aq_ldiv_B!(X, lu, B, UMFPACK_A)
960-
ldiv!(X::StridedVecOrMat{Tb}, translu::TransType{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{Tb}) where {Tb<:Complex} =
960+
ldiv!(X::StridedVecOrMat{Tb}, translu::TransposeFact{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{Tb}) where {Tb<:Complex} =
961961
(lu = translu.parent; _Aq_ldiv_B!(X, lu, B, UMFPACK_Aat))
962-
ldiv!(X::StridedVecOrMat{Tb}, adjlu::AdjType{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{Tb}) where {Tb<:Complex} =
962+
ldiv!(X::StridedVecOrMat{Tb}, adjlu::AdjointFact{Float64,<:UmfpackLU{Float64}}, B::StridedVecOrMat{Tb}) where {Tb<:Complex} =
963963
(lu = adjlu.parent; _Aq_ldiv_B!(X, lu, B, UMFPACK_At))
964964

965965
function _Aq_ldiv_B!(X::StridedVecOrMat, lu::UmfpackLU, B::StridedVecOrMat, transposeoptype)

‎test/umfpack.jl

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,7 @@ using Random
77
using SparseArrays
88
using Serialization
99
using LinearAlgebra:
10-
I, det, issuccess, ldiv!, lu, lu!, Adjoint, Transpose, SingularException, Diagonal, logabsdet
10+
LinearAlgebra, I, det, issuccess, ldiv!, lu, lu!, Transpose, SingularException, Diagonal, logabsdet
1111
using SparseArrays: nnz, sparse, sprand, sprandn, SparseMatrixCSC, UMFPACK, increment!
1212
if Base.USE_GPL_LIBS
1313
function umfpack_report(l::UMFPACK.UmfpackLU)
@@ -16,6 +16,10 @@ function umfpack_report(l::UMFPACK.UmfpackLU)
1616
return
1717
end
1818

19+
const TransposeFact = isdefined(LinearAlgebra, :TransposeFactorization) ?
20+
LinearAlgebra.TransposeFactorization :
21+
Transpose
22+
1923
for itype in UMFPACK.UmfpackIndexTypes
2024
sol_r = Symbol(UMFPACK.umf_nm("solve", :Float64, itype))
2125
sol_c = Symbol(UMFPACK.umf_nm("solve", :ComplexF64, itype))
@@ -218,6 +222,7 @@ end
218222
@test y ≈ x
219223

220224
@test A'*x ≈ b
225+
@test transpose(lua) isa TransposeFact
221226
x = transpose(lua) \ b
222227
@test x ≈ float([1:5;])
223228

0 commit comments

Comments
 (0)