Adjust to AdjointFactorization (#313)

JuliaLinearAlgebra · Feb 12, 2023 · 6c66f3e · 6c66f3e · dkarrasch · Feb 12, 2023
1 parent 992ab32
commit 6c66f3e
Show file tree

Hide file tree

Showing 5 changed files with 27 additions and 13 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,6 +1,6 @@
 name = "BandedMatrices"
 uuid = "aae01518-5342-5314-be14-df237901396f"
-version = "0.17.11"
+version = "0.17.12"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

diff --git a/src/BandedMatrices.jl b/src/BandedMatrices.jl
@@ -40,6 +40,12 @@ import FillArrays: AbstractFill, getindex_value, _broadcasted_zeros, unique_valu
 
 const libblas = Base.libblas_name
 const liblapack = Base.liblapack_name
+const AdjointFact = isdefined(LinearAlgebra, :AdjointFactorization) ?
+    LinearAlgebra.AdjointFactorization :
+    Adjoint
+const TransposeFact = isdefined(LinearAlgebra, :TransposeFactorization) ?
+    LinearAlgebra.TransposeFactorization :
+    Transpose
 
 export BandedMatrix,
        bandrange,

diff --git a/src/banded/BandedLU.jl b/src/banded/BandedLU.jl
@@ -32,8 +32,10 @@ Base.iterate(S::BandedLU, ::Val{:U}) = (S.U, Val(:p))
 Base.iterate(S::BandedLU, ::Val{:p}) = (S.p, Val(:done))
 Base.iterate(S::BandedLU, ::Val{:done}) = nothing
 
-adjoint(F::BandedLU) = Adjoint(F)
-transpose(F::BandedLU) = Transpose(F)
+if !(isdefined(LinearAlgebra, :AdjointFactorization)) # VERSION < v"1.10-"
+    adjoint(F::BandedLU) = Adjoint(F)
+end
+transpose(F::BandedLU) = TransposeFact(F)
 
 lu(S::BandedLU) = S
 

diff --git a/src/banded/linalg.jl b/src/banded/linalg.jl
@@ -38,38 +38,38 @@ function ldiv!(A::BandedLU{T}, B::AbstractVecOrMat{Complex{T}}) where T<:Real
     B .= a .+ im.*b
 end
 
-function ldiv!(transA::Transpose{T,<:BandedLU{T,<:BandedMatrix}}, B::StridedVecOrMat{T}) where {T<:BlasFloat}
+function ldiv!(transA::TransposeFact{T,<:BandedLU{T,<:BandedMatrix}}, B::StridedVecOrMat{T}) where {T<:BlasFloat}
     A = transA.parent
     m = size(A.factors,1)
     l,u = bandwidths(A.factors)
     data = bandeddata(A.factors)
     LAPACK.gbtrs!('T', l, u-l, m, data, A.ipiv, B)
 end
 
-function ldiv!(transA::Transpose{<:Any,<:BandedLU{<:Any,<:BandedMatrix}}, B::AbstractVecOrMat)
+function ldiv!(transA::TransposeFact{<:Any,<:BandedLU{<:Any,<:BandedMatrix}}, B::AbstractVecOrMat)
     A = transA.parent
     ldiv!(transpose(UnitLowerTriangular(A.factors)), ldiv!(transpose(UpperTriangular(A.factors)), B))
     _apply_inverse_ipiv_rows!(A, B)
 end
 
-ldiv!(adjF::Adjoint{T,<:BandedLU{T,<:BandedMatrix}}, B::AbstractVecOrMat{T}) where {T<:Real} =
+ldiv!(adjF::AdjointFact{T,<:BandedLU{T,<:BandedMatrix}}, B::AbstractVecOrMat{T}) where {T<:Real} =
     (F = adjF.parent; ldiv!(transpose(F), B))
-function ldiv!(adjA::Adjoint{T,<:BandedLU{T,<:BandedMatrix}}, B::StridedVecOrMat{T}) where {T<:BlasComplex}
+function ldiv!(adjA::AdjointFact{T,<:BandedLU{T,<:BandedMatrix}}, B::StridedVecOrMat{T}) where {T<:BlasComplex}
     A = adjA.parent
     m = size(A.factors,1)
     l,u = bandwidths(A.factors)
     data = bandeddata(A.factors)
     LAPACK.gbtrs!('C', l, u-l, m, data, A.ipiv, B)
 end
 
-function ldiv!(adjA::Adjoint{<:Any,<:BandedLU{<:Any,<:BandedMatrix}}, B::AbstractVecOrMat)
+function ldiv!(adjA::AdjointFact{<:Any,<:BandedLU{<:Any,<:BandedMatrix}}, B::AbstractVecOrMat)
     error("Implement")
     A = adjA.parent
     ldiv!(adjoint(UnitLowerTriangular(A.factors)), ldiv!(adjoint(UpperTriangular(A.factors)), B))
     _apply_inverse_ipiv!(A, B)
 end
 
-\(A::Adjoint{<:Any,<:BandedLU}, B::Adjoint{<:Any,<:AbstractVecOrMat}) = A \ copy(B)
-\(A::Transpose{<:Any,<:BandedLU}, B::Transpose{<:Any,<:AbstractVecOrMat}) = A \ copy(B)
+\(A::AdjointFact{<:Any,<:BandedLU}, B::Adjoint{<:Any,<:AbstractVecOrMat}) = A \ copy(B)
+\(A::TransposeFact{<:Any,<:BandedLU}, B::Transpose{<:Any,<:AbstractVecOrMat}) = A \ copy(B)
 
 _factorize(::AbstractBandedLayout, _, A) = size(A,1) == size(A,2) ? lu(A) : qr(A)
diff --git a/src/symbanded/SplitCholesky.jl b/src/symbanded/SplitCholesky.jl
@@ -45,7 +45,13 @@ function splitcholesky!(::SymmetricLayout{<:BandedColumnMajor},
     SplitCholesky(A, A.uplo)
 end
 
-adjoint(S::SplitCholesky) = Adjoint(S)
+if !(isdefined(LinearAlgebra, :AdjointFactorization)) # VERSION < v"1.10-"
+    adjoint(S::SplitCholesky) = Adjoint(S)
+else
+    transpose(S::SplitCholesky{<:Real}) = S'
+    transpose(::SplitCholesky) =
+        throw(ArgumentError("transpose of SplitCholesky decomposition is not supported, consider using adjoint"))
+end
 
 function lmul!(S::SplitCholesky{T,Symmetric{T,M}}, B::AbstractVecOrMat{T}) where {T,M<:BandedMatrix{T}}
     require_one_based_indexing(B)
@@ -103,7 +109,7 @@ function ldiv!(S::SplitCholesky{T,Symmetric{T,M}}, B::AbstractVecOrMat{T}) where
     B
 end
 
-function lmul!(S::Adjoint{T,SplitCholesky{T,Symmetric{T,M}}}, B::AbstractVecOrMat{T}) where {T,M<:BandedMatrix{T}}
+function lmul!(S::AdjointFact{T,SplitCholesky{T,Symmetric{T,M}}}, B::AbstractVecOrMat{T}) where {T,M<:BandedMatrix{T}}
     require_one_based_indexing(B)
     n, nrhs = size(B, 1), size(B, 2)
     if size(S, 1) != n
@@ -138,7 +144,7 @@ function lmul!(S::Adjoint{T,SplitCholesky{T,Symmetric{T,M}}}, B::AbstractVecOrMa
     B
 end
 
-function ldiv!(S::Adjoint{T,SplitCholesky{T,Symmetric{T,M}}}, B::AbstractVecOrMat{T}) where {T,M<:BandedMatrix{T}}
+function ldiv!(S::AdjointFact{T,SplitCholesky{T,Symmetric{T,M}}}, B::AbstractVecOrMat{T}) where {T,M<:BandedMatrix{T}}
     require_one_based_indexing(B)
     n, nrhs = size(B, 1), size(B, 2)
     if size(S, 1) != n