From 8112840502344da874a920dc37d83601de0e60e3 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Wed, 3 Aug 2022 21:22:32 +0200 Subject: [PATCH 01/16] Adjustment to `AdjointQ` in LinearAlgebra --- src/ArrayLayouts.jl | 7 ++++++ src/factorizations.jl | 50 +++++++++++++++++++++---------------------- src/lmul.jl | 4 ++-- src/memorylayout.jl | 3 +++ src/mul.jl | 4 ++++ 5 files changed, 41 insertions(+), 27 deletions(-) diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index 9d3c2c4..5cc8cf4 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -39,6 +39,13 @@ import LinearAlgebra: AbstractTriangular, AbstractQ, QRCompactWYQ, QRPackedQ, ch norm2, norm1, normInf, normMinusInf, qr, lu, qr!, lu!, AdjOrTrans, HermOrSym, AdjointAbsVec, TransposeAbsVec, cholcopy, checknonsingular, _apply_ipiv_rows!, ipiv2perm, RealHermSymComplexHerm, chkfullrank +if VERSION >= v"1.9-" + using LinearAlgebra: AdjointQ + AdjointQtype{T} = AdjointQ{T} +else + AdjointQtype{T} = Adjoint{T,<:AbstractQ} +end + import LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex import FillArrays: AbstractFill, getindex_value, axes_print_matrix_row, _copy_oftype diff --git a/src/factorizations.jl b/src/factorizations.jl index 4426887..1c27739 100644 --- a/src/factorizations.jl +++ b/src/factorizations.jl @@ -123,7 +123,7 @@ adjointlayout(::Type, ::QRPackedQLayout{SLAY,TLAY}) where {SLAY,TLAY} = AdjQRPac adjointlayout(::Type, ::QRCompactWYQLayout{SLAY,TLAY}) where {SLAY,TLAY} = AdjQRCompactWYQLayout{SLAY,TLAY}() copy(M::Lmul{<:AbstractQLayout}) = copyto!(similar(M), M) -mulreduce(M::Mul{<:AbstractQLayout,<:AbstractQLayout}) = MulAdd(M) +mulreduce(M::Mul{<:AbstractQLayout,<:AbstractQLayout}) = Lmul(M) mulreduce(M::Mul{<:AbstractQLayout}) = Lmul(M) mulreduce(M::Mul{<:Any,<:AbstractQLayout}) = Rmul(M) mulreduce(M::Mul{<:TriangularLayout,<:AbstractQLayout}) = Rmul(M) @@ -151,17 +151,17 @@ materialize!(M::Rmul{LAY}) where LAY<:AbstractQLayout = error("Overload material materialize!(M::Ldiv{<:AbstractQLayout}) = materialize!(Lmul(M.A',M.B)) -materialize!(M::Lmul{<:QRPackedQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractColumnMajor,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = +materialize!(M::Lmul{<:QRPackedQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractColumnMajor,<:AbstractQ{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = LAPACK.ormqr!('L','N',M.A.factors,M.A.τ,M.B) -materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractColumnMajor,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = +materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractColumnMajor,<:AbstractQ{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = LAPACK.gemqrt!('L','N',M.A.factors,M.A.T,M.B) -function materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractRowMajor,<:AbstractMatrix{T},<:AbstractMatrix{T}}) where T<:BlasReal +function materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:AbstractRowMajor,<:AbstractQ{T},<:AbstractMatrix{T}}) where T<:BlasReal LAPACK.gemqrt!('R','T',M.A.factors,M.A.T,transpose(M.B)) M.B end -function materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:ConjLayout{<:AbstractRowMajor},<:AbstractMatrix{T},<:AbstractMatrix{T}}) where T<:BlasComplex +function materialize!(M::Lmul{<:QRCompactWYQLayout{<:AbstractColumnMajor,<:AbstractColumnMajor},<:ConjLayout{<:AbstractRowMajor},<:AbstractQ{T},<:AbstractMatrix{T}}) where T<:BlasComplex LAPACK.gemqrt!('R','C',M.A.factors,M.A.T,(M.B)') M.B end @@ -196,18 +196,18 @@ end ### QcB -materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = - (A = M.A.parent; LAPACK.ormqr!('L','T',A.factors,A.τ,M.B)) -materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = - (A = M.A.parent; LAPACK.ormqr!('L','C',A.factors,A.τ,M.B)) -materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = - (A = M.A.parent; LAPACK.gemqrt!('L','T',A.factors,A.T,M.B)) -materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = - (A = M.A.parent; LAPACK.gemqrt!('L','C',A.factors,A.T,M.B)) +materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = + (A = parent(M.A); LAPACK.ormqr!('L','T',A.factors,A.τ,M.B)) +materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = + (A = parent(M.A); LAPACK.ormqr!('L','C',A.factors,A.τ,M.B)) +materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = + (A = parent(M.A); LAPACK.gemqrt!('L','T',A.factors,A.T,M.B)) +materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = + (A = parent(M.A); LAPACK.gemqrt!('L','C',A.factors,A.T,M.B)) function materialize!(M::Lmul{<:AdjQRPackedQLayout}) adjA,B = M.A, M.B require_one_based_indexing(B) - A = adjA.parent + A = parent(adjA) mA, nA = size(A.factors) mB, nB = size(B,1), size(B,2) if mA != mB @@ -233,9 +233,9 @@ function materialize!(M::Lmul{<:AdjQRPackedQLayout}) end ## AQ -materialize!(M::Rmul{<:AbstractStridedLayout,<:QRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasFloat = +materialize!(M::Rmul{<:AbstractStridedLayout,<:QRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractQ{T}}) where T<:BlasFloat = LAPACK.ormqr!('R', 'N', M.B.factors, M.B.τ, M.A) -materialize!(M::Rmul{<:AbstractStridedLayout,<:QRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasFloat = +materialize!(M::Rmul{<:AbstractStridedLayout,<:QRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractQ{T}}) where T<:BlasFloat = LAPACK.gemqrt!('R','N', M.B.factors, M.B.T, M.A) function materialize!(M::Rmul{<:Any,<:QRPackedQLayout}) A,Q = M.A,M.B @@ -264,17 +264,17 @@ function materialize!(M::Rmul{<:Any,<:QRPackedQLayout}) end ### AQc -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasReal = - (B = M.B.parent; LAPACK.ormqr!('R','T',B.factors,B.τ,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasComplex = - (B = M.B.parent; LAPACK.ormqr!('R','C',B.factors,B.τ,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasReal = - (B = M.B.parent; LAPACK.gemqrt!('R','T',B.factors,B.T,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasComplex = - (B = M.B.parent; LAPACK.gemqrt!('R','C',B.factors,B.T,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasReal = + (B = parent(M.B); LAPACK.ormqr!('R','T',B.factors,B.τ,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasComplex = + (B = parent(M.B); LAPACK.ormqr!('R','C',B.factors,B.τ,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasReal = + (B = parent(M.B); LAPACK.gemqrt!('R','T',B.factors,B.T,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasComplex = + (B = parent(M.B); LAPACK.gemqrt!('R','C',B.factors,B.T,M.A)) function materialize!(M::Rmul{<:Any,<:AdjQRPackedQLayout}) A,adjQ = M.A,M.B - Q = adjQ.parent + Q = parent(adjQ) mQ, nQ = size(Q.factors) mA, nA = size(A,1), size(A,2) if nA != mQ diff --git a/src/lmul.jl b/src/lmul.jl index 0143b67..253062d 100644 --- a/src/lmul.jl +++ b/src/lmul.jl @@ -27,10 +27,10 @@ end -const MatLmulVec{StyleA,StyleB} = Lmul{StyleA,StyleB,<:AbstractMatrix,<:AbstractVector} +const MatLmulVec{StyleA,StyleB} = Lmul{StyleA,StyleB,<:Union{AbstractMatrix,AbstractQ},<:AbstractVector} const MatLmulMat{StyleA,StyleB} = Lmul{StyleA,StyleB,<:AbstractMatrix,<:AbstractMatrix} -const BlasMatLmulVec{StyleA,StyleB,T<:BlasFloat} = Lmul{StyleA,StyleB,<:AbstractMatrix{T},<:AbstractVector{T}} +const BlasMatLmulVec{StyleA,StyleB,T<:BlasFloat} = Lmul{StyleA,StyleB,<:Union{AbstractMatrix{T},AbstractQ{T}},<:AbstractVector{T}} const BlasMatLmulMat{StyleA,StyleB,T<:BlasFloat} = Lmul{StyleA,StyleB,<:AbstractMatrix{T},<:AbstractMatrix{T}} const MatRmulMat{StyleA,StyleB} = Rmul{StyleA,StyleB,<:AbstractMatrix,<:AbstractMatrix} diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 79f8290..348f5aa 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -281,6 +281,9 @@ struct DualLayout{ML<:MemoryLayout} <: MemoryLayout end MemoryLayout(::Type{Transpose{T,P}}) where {T,P} = transposelayout(MemoryLayout(P)) MemoryLayout(::Type{Adjoint{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) +if VERSION >= v"1.9-" + MemoryLayout(::Type{AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) +end MemoryLayout(::Type{AdjointAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(adjointlayout(T,MemoryLayout(P)))}() MemoryLayout(::Type{TransposeAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(transposelayout(MemoryLayout(P)))}() diff --git a/src/mul.jl b/src/mul.jl index be0ab69..3ecd951 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -196,6 +196,10 @@ macro layoutmul(Typ) Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) + if VERSION >= v"1.9-" + Base.:*(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) + Base.:*(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) + end end for Struc in (:AbstractTriangular, :Diagonal) ret = quote From 8d349edd59052f5d7f5e6a3d078239c5468a4704 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 24 Sep 2022 13:47:59 +0200 Subject: [PATCH 02/16] Update ArrayLayouts.jl --- src/ArrayLayouts.jl | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index 5cc8cf4..be5dd14 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -39,12 +39,7 @@ import LinearAlgebra: AbstractTriangular, AbstractQ, QRCompactWYQ, QRPackedQ, ch norm2, norm1, normInf, normMinusInf, qr, lu, qr!, lu!, AdjOrTrans, HermOrSym, AdjointAbsVec, TransposeAbsVec, cholcopy, checknonsingular, _apply_ipiv_rows!, ipiv2perm, RealHermSymComplexHerm, chkfullrank -if VERSION >= v"1.9-" - using LinearAlgebra: AdjointQ - AdjointQtype{T} = AdjointQ{T} -else - AdjointQtype{T} = Adjoint{T,<:AbstractQ} -end +AdjointQtype{T} = isdefined(LinearAlgebra, :AdjointQ) ? LinearAlgebra.AdjointQ{T} : Adjoint{T,<:AbstractQ} import LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex From a7f2610e3ccebe436774cbc4416a26760de0b3a4 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 24 Sep 2022 13:50:16 +0200 Subject: [PATCH 03/16] Update memorylayout.jl --- src/memorylayout.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 348f5aa..242b5a7 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -281,7 +281,7 @@ struct DualLayout{ML<:MemoryLayout} <: MemoryLayout end MemoryLayout(::Type{Transpose{T,P}}) where {T,P} = transposelayout(MemoryLayout(P)) MemoryLayout(::Type{Adjoint{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) -if VERSION >= v"1.9-" +if isdefined(LinearAlgebra, :AdjointQ) MemoryLayout(::Type{AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) end MemoryLayout(::Type{AdjointAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(adjointlayout(T,MemoryLayout(P)))}() From 64d9730b2f3f2bb7d76930f7d4b4c8e87e9ab423 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 24 Sep 2022 13:51:06 +0200 Subject: [PATCH 04/16] Update mul.jl --- src/mul.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/mul.jl b/src/mul.jl index 3ecd951..d157056 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -196,7 +196,7 @@ macro layoutmul(Typ) Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) - if VERSION >= v"1.9-" + if isdefined(LinearAlgebra, AdjointQ) Base.:*(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) end From 20339fc876e4e61ac0a782ed7725fd64c536b073 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 24 Sep 2022 13:55:54 +0200 Subject: [PATCH 05/16] Update mul.jl --- src/mul.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/mul.jl b/src/mul.jl index d157056..f0ad087 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -196,7 +196,7 @@ macro layoutmul(Typ) Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) - if isdefined(LinearAlgebra, AdjointQ) + if isdefined(LinearAlgebra, :AdjointQ) Base.:*(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) end From 94ef47dbde4b8f73af424408254186418e00e20b Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 24 Sep 2022 15:08:53 +0200 Subject: [PATCH 06/16] bump patch number --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 18fa539..39c8164 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ArrayLayouts" uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" authors = ["Sheehan Olver "] -version = "0.8.11" +version = "0.8.12" [deps] FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" From 941f8c0634dd44052555346475ad2d7df2ff23e4 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Tue, 4 Oct 2022 19:15:59 +0200 Subject: [PATCH 07/16] further adjustments --- src/memorylayout.jl | 2 +- src/mul.jl | 5 +++++ 2 files changed, 6 insertions(+), 1 deletion(-) diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 242b5a7..6cfd72b 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -282,7 +282,7 @@ struct DualLayout{ML<:MemoryLayout} <: MemoryLayout end MemoryLayout(::Type{Transpose{T,P}}) where {T,P} = transposelayout(MemoryLayout(P)) MemoryLayout(::Type{Adjoint{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) if isdefined(LinearAlgebra, :AdjointQ) - MemoryLayout(::Type{AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) + MemoryLayout(::Type{LinearAlgebra.AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) end MemoryLayout(::Type{AdjointAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(adjointlayout(T,MemoryLayout(P)))}() MemoryLayout(::Type{TransposeAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(transposelayout(MemoryLayout(P)))}() diff --git a/src/mul.jl b/src/mul.jl index f0ad087..dfbfc3b 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -196,6 +196,11 @@ macro layoutmul(Typ) Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) + # disambiguation + Base.:*(A::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}, B::$Typ) = + ArrayLayouts.mul(A,B) + Base.:*(A::$Typ, B::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}) = + ArrayLayouts.mul(A,B) if isdefined(LinearAlgebra, :AdjointQ) Base.:*(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) Base.:*(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) From a9c6e5593afe411e7bdbe1699e6c2e8b9ed50700 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Wed, 2 Nov 2022 10:53:33 +0100 Subject: [PATCH 08/16] mild clean-up --- Project.toml | 3 +- src/ArrayLayouts.jl | 69 +++++++++++++++------------------------- src/ldiv.jl | 2 ++ test/test_layoutarray.jl | 2 +- test/test_muladd.jl | 2 +- 5 files changed, 31 insertions(+), 47 deletions(-) diff --git a/Project.toml b/Project.toml index 39c8164..3c0fb86 100644 --- a/Project.toml +++ b/Project.toml @@ -14,10 +14,9 @@ julia = "1.6" [extras] Base64 = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" -Compat = "34da2185-b29b-5c13-b0c7-acf172513d20" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Base64", "Compat", "Random", "StableRNGs", "Test"] +test = ["Base64", "Random", "StableRNGs", "Test"] diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index be5dd14..3982b7d 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -1,59 +1,41 @@ module ArrayLayouts using Base: _typed_hcat using Base, Base.Broadcast, LinearAlgebra, FillArrays, SparseArrays -import LinearAlgebra.BLAS - -import Base: AbstractArray, AbstractMatrix, AbstractVector, - ReinterpretArray, ReshapedArray, AbstractCartesianIndex, Slice, - RangeIndex, BroadcastStyle, copyto!, length, broadcastable, axes, - getindex, eltype, tail, IndexStyle, IndexLinear, getproperty, - *, +, -, /, \, ==, isinf, isfinite, sign, angle, show, isless, - fld, cld, div, min, max, minimum, maximum, mod, - <, ≤, >, ≥, promote_rule, convert, copy, - size, step, isempty, length, first, last, ndims, - getindex, setindex!, intersect, @_inline_meta, inv, - sort, sort!, issorted, sortperm, diff, cumsum, sum, in, broadcast, - eltype, parent, real, imag, - conj, transpose, adjoint, permutedims, vec, - exp, log, sqrt, cos, sin, tan, csc, sec, cot, - cosh, sinh, tanh, csch, sech, coth, - acos, asin, atan, acsc, asec, acot, - acosh, asinh, atanh, acsch, asech, acoth, (:), - AbstractMatrix, AbstractArray, checkindex, unsafe_length, OneTo, one, zero, - to_shape, _sub2ind, print_matrix, print_matrix_row, print_matrix_vdots, - checkindex, Slice, @propagate_inbounds, @_propagate_inbounds_meta, - _in_range, _range, Ordered, - ArithmeticWraps, floatrange, reverse, unitrange_last, - AbstractArray, AbstractVector, axes, (:), _sub2ind_recurse, broadcast, promote_eltypeof, - similar, @_gc_preserve_end, @_gc_preserve_begin, - @nexprs, @ncall, @ntuple, tuple_type_tail, - all, any, isbitsunion, issubset, replace_in_print_matrix, replace_with_centered_mark, - strides, unsafe_convert, first_index, unalias, union - -import Base.Broadcast: BroadcastStyle, AbstractArrayStyle, Broadcasted, broadcasted, - combine_eltypes, DefaultArrayStyle, instantiate, materialize, - materialize!, eltypes - -import LinearAlgebra: AbstractTriangular, AbstractQ, QRCompactWYQ, QRPackedQ, checksquare, pinv, - fill!, tilebufsize, factorize, qr, lu, cholesky, - norm2, norm1, normInf, normMinusInf, qr, lu, qr!, lu!, AdjOrTrans, HermOrSym, AdjointAbsVec, - TransposeAbsVec, cholcopy, checknonsingular, _apply_ipiv_rows!, ipiv2perm, RealHermSymComplexHerm, chkfullrank +using LinearAlgebra.BLAS + +using Base: AbstractCartesianIndex, OneTo, RangeIndex, ReinterpretArray, ReshapedArray, + Slice, tuple_type_tail, unalias, + @propagate_inbounds, @_propagate_inbounds_meta + +import Base: axes, size, length, eltype, ndims, first, last, diff, isempty, union, sort!, + ==, *, +, -, /, \, copy, copyto!, similar, getindex, strides, + unsafe_convert + +using Base.Broadcast: Broadcasted + +import Base.Broadcast: BroadcastStyle, broadcastable, instantiate, materialize, materialize! + +using LinearAlgebra: AbstractTriangular, AbstractQ, QRCompactWYQ, QRPackedQ, checksquare, + pinv, tilebufsize, cholcopy, + norm2, norm1, normInf, normMinusInf, + AdjOrTrans, HermOrSym, RealHermSymComplexHerm, AdjointAbsVec, TransposeAbsVec, + checknonsingular, _apply_ipiv_rows!, ipiv2perm, chkfullrank AdjointQtype{T} = isdefined(LinearAlgebra, :AdjointQ) ? LinearAlgebra.AdjointQ{T} : Adjoint{T,<:AbstractQ} -import LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex +using LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex -import FillArrays: AbstractFill, getindex_value, axes_print_matrix_row, _copy_oftype +using FillArrays: AbstractFill, getindex_value, axes_print_matrix_row, _copy_oftype -import Base: require_one_based_indexing +using Base: require_one_based_indexing export materialize, materialize!, MulAdd, muladd!, Ldiv, Rdiv, Lmul, Rmul, Dot, - lmul, rmul, mul, ldiv, rdiv, mul, MemoryLayout, AbstractStridedLayout, + lmul, mul, ldiv, rdiv, mul, MemoryLayout, AbstractStridedLayout, DenseColumnMajor, ColumnMajor, ZerosLayout, FillLayout, AbstractColumnMajor, RowMajor, AbstractRowMajor, UnitStride, DiagonalLayout, ScalarLayout, SymTridiagonalLayout, TridiagonalLayout, BidiagonalLayout, HermitianLayout, SymmetricLayout, TriangularLayout, UnknownLayout, AbstractBandedLayout, ApplyBroadcastStyle, ConjLayout, AbstractFillLayout, DualLayout, - colsupport, rowsupport, layout_getindex, QLayout, LayoutArray, LayoutMatrix, LayoutVector, + colsupport, rowsupport, layout_getindex, AbstractQLayout, LayoutArray, LayoutMatrix, LayoutVector, RangeCumsum if VERSION < v"1.7-" @@ -83,7 +65,8 @@ abstract type LayoutArray{T,N} <: AbstractArray{T,N} end const LayoutMatrix{T} = LayoutArray{T,2} const LayoutVector{T} = LayoutArray{T,1} -## TODO: Following are type piracy whch may be removed in Julia v1.5 +## TODO: Following are type piracy which may be removed in Julia v1.5 +## No, it can't, because strides(::AdjointAbsMat) is defined only for real eltype! _transpose_strides(a) = (a,1) _transpose_strides(a,b) = (b,a) strides(A::Adjoint) = _transpose_strides(strides(parent(A))...) diff --git a/src/ldiv.jl b/src/ldiv.jl index f03cd71..a748e41 100644 --- a/src/ldiv.jl +++ b/src/ldiv.jl @@ -133,6 +133,8 @@ macro _layoutldiv(Typ) LinearAlgebra.ldiv!(A::$Typ, x::StridedMatrix) = ArrayLayouts.ldiv!(A,x) LinearAlgebra.ldiv!(A::Factorization, x::$Typ) = ArrayLayouts.ldiv!(A,x) + LinearAlgebra.ldiv!(A::LU, x::$Typ) = ArrayLayouts.ldiv!(A,x) + LinearAlgebra.ldiv!(A::Cholesky, x::$Typ) = ArrayLayouts.ldiv!(A,x) LinearAlgebra.ldiv!(A::Bidiagonal, B::$Typ) = ArrayLayouts.ldiv!(A,B) diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index a7008bd..dd78487 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -1,5 +1,5 @@ using ArrayLayouts, LinearAlgebra, FillArrays, Base64, Test -import ArrayLayouts: sub_materialize, MemoryLayout, ColumnNorm, RowMaximum, CRowMaximum +using ArrayLayouts: sub_materialize, MemoryLayout, ColumnNorm, RowMaximum, CRowMaximum struct MyMatrix <: LayoutMatrix{Float64} diff --git a/test/test_muladd.jl b/test/test_muladd.jl index b586a82..731a990 100644 --- a/test/test_muladd.jl +++ b/test/test_muladd.jl @@ -1,5 +1,5 @@ using ArrayLayouts, FillArrays, Random, StableRNGs, LinearAlgebra, Test -import ArrayLayouts: DenseColumnMajor, AbstractStridedLayout, AbstractColumnMajor, DiagonalLayout, mul, Mul, zero! +using ArrayLayouts: DenseColumnMajor, AbstractStridedLayout, AbstractColumnMajor, DiagonalLayout, mul, Mul, zero! Random.seed!(0) @testset "Multiplication" begin From 00cfa2d29bb5835b0d47b969de844d46720201e8 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Wed, 2 Nov 2022 11:42:48 +0100 Subject: [PATCH 09/16] import getproperty --- src/ArrayLayouts.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index 3982b7d..044115c 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -8,7 +8,7 @@ using Base: AbstractCartesianIndex, OneTo, RangeIndex, ReinterpretArray, Reshape @propagate_inbounds, @_propagate_inbounds_meta import Base: axes, size, length, eltype, ndims, first, last, diff, isempty, union, sort!, - ==, *, +, -, /, \, copy, copyto!, similar, getindex, strides, + ==, *, +, -, /, \, copy, copyto!, similar, getproperty, getindex, strides, unsafe_convert using Base.Broadcast: Broadcasted From 1cd53244a88a866032fe7642298d4f2b87a41ded Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Wed, 2 Nov 2022 12:21:02 +0100 Subject: [PATCH 10/16] rm specifiers of imported/exported functions --- src/ArrayLayouts.jl | 28 +++++------ src/cumsum.jl | 4 +- src/factorizations.jl | 12 ++--- src/ldiv.jl | 46 +++++++++--------- src/memorylayout.jl | 18 +++---- src/mul.jl | 110 +++++++++++++++++++++--------------------- src/muladd.jl | 2 +- 7 files changed, 110 insertions(+), 110 deletions(-) diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index 044115c..d8117c5 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -9,7 +9,7 @@ using Base: AbstractCartesianIndex, OneTo, RangeIndex, ReinterpretArray, Reshape import Base: axes, size, length, eltype, ndims, first, last, diff, isempty, union, sort!, ==, *, +, -, /, \, copy, copyto!, similar, getproperty, getindex, strides, - unsafe_convert + reverse, unsafe_convert using Base.Broadcast: Broadcasted @@ -130,16 +130,16 @@ Base.@propagate_inbounds layout_getindex(A::AbstractArray, I::CartesianIndex) = macro _layoutgetindex(Typ) esc(quote - @inline Base.getindex(A::$Typ, kr::Colon, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::Colon, jr::AbstractUnitRange) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::AbstractUnitRange, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::AbstractUnitRange, jr::AbstractUnitRange) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::AbstractVector, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::Colon, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::Colon, jr::Integer) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::AbstractVector, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::Integer, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) - @inline Base.getindex(A::$Typ, kr::Integer, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Colon, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Colon, jr::AbstractUnitRange) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::AbstractUnitRange, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::AbstractUnitRange, jr::AbstractUnitRange) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::AbstractVector, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Colon, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Colon, jr::Integer) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::AbstractVector, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Integer, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr) + @inline getindex(A::$Typ, kr::Integer, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr) end) end @@ -228,7 +228,7 @@ end _copyto!(_, _, dest::AbstractArray{T,N}, src::AbstractArray{V,N}) where {T,V,N} = - Base.invoke(copyto!, Tuple{AbstractArray{T,N},AbstractArray{V,N}}, dest, src) + invoke(copyto!, Tuple{AbstractArray{T,N},AbstractArray{V,N}}, dest, src) _copyto!(dest, src) = _copyto!(MemoryLayout(dest), MemoryLayout(src), dest, src) @@ -266,7 +266,7 @@ function zero!(_, A::AbstractArray{<:AbstractArray}) A end -_norm(_, A, p) = Base.invoke(norm, Tuple{Any,Real}, A, p) +_norm(_, A, p) = invoke(norm, Tuple{Any,Real}, A, p) LinearAlgebra.norm(A::LayoutArray, p::Real=2) = _norm(MemoryLayout(A), A, p) LinearAlgebra.norm(A::SubArray{<:Any,N,<:LayoutArray}, p::Real=2) where N = _norm(MemoryLayout(A), A, p) @@ -338,7 +338,7 @@ layout_replace_in_print_matrix(_, A, i, j, s) = i in colsupport(A,j) ? s : Base.replace_with_centered_mark(s) Base.replace_in_print_matrix(A::Union{LayoutVector, - LayoutMatrix, + LayoutMatrix, UpperTriangular{<:Any,<:LayoutMatrix}, UnitUpperTriangular{<:Any,<:LayoutMatrix}, LowerTriangular{<:Any,<:LayoutMatrix}, diff --git a/src/cumsum.jl b/src/cumsum.jl index 5de24fa..9f225cf 100644 --- a/src/cumsum.jl +++ b/src/cumsum.jl @@ -31,7 +31,7 @@ last(r::RangeCumsum) = sum(r.range) diff(r::RangeCumsum) = r.range[2:end] isempty(r::RangeCumsum) = isempty(r.range) -union(a::RangeCumsum{<:Any,<:Base.OneTo}, b::RangeCumsum{<:Any,<:Base.OneTo}) = - RangeCumsum(Base.OneTo(max(last(a.range),last(b.range)))) +union(a::RangeCumsum{<:Any,<:OneTo}, b::RangeCumsum{<:Any,<:OneTo}) = + RangeCumsum(OneTo(max(last(a.range), last(b.range)))) sort!(a::RangeCumsum{<:Any,<:AbstractUnitRange}) = a \ No newline at end of file diff --git a/src/factorizations.jl b/src/factorizations.jl index 1c27739..26608fb 100644 --- a/src/factorizations.jl +++ b/src/factorizations.jl @@ -300,12 +300,12 @@ function materialize!(M::Rmul{<:Any,<:AdjQRPackedQLayout}) end -__qr(layout, lengths, A; kwds...) = Base.invoke(qr, Tuple{AbstractMatrix{eltype(A)}}, A; kwds...) +__qr(layout, lengths, A; kwds...) = invoke(qr, Tuple{AbstractMatrix{eltype(A)}}, A; kwds...) _qr(layout, axes, A; kwds...) = __qr(layout, map(length, axes), A; kwds...) -_qr(layout, axes, A, pivot::P; kwds...) where P = Base.invoke(qr, Tuple{AbstractMatrix{eltype(A)},P}, A, pivot; kwds...) +_qr(layout, axes, A, pivot::P; kwds...) where P = invoke(qr, Tuple{AbstractMatrix{eltype(A)},P}, A, pivot; kwds...) _qr!(layout, axes, A, args...; kwds...) = error("Overload _qr!(::$(typeof(layout)), axes, A)") -_lu(layout, axes, A; kwds...) = Base.invoke(lu, Tuple{AbstractMatrix{eltype(A)}}, A; kwds...) -_lu(layout, axes, A, pivot::P; kwds...) where P = Base.invoke(lu, Tuple{AbstractMatrix{eltype(A)},P}, A, pivot; kwds...) +_lu(layout, axes, A; kwds...) = invoke(lu, Tuple{AbstractMatrix{eltype(A)}}, A; kwds...) +_lu(layout, axes, A, pivot::P; kwds...) where P = invoke(lu, Tuple{AbstractMatrix{eltype(A)},P}, A, pivot; kwds...) _lu!(layout, axes, A, args...; kwds...) = error("Overload _lu!(::$(typeof(layout)), axes, A)") _cholesky(layout, axes, A, ::CNoPivot=CNoPivot(); check::Bool = true) = cholesky!(cholcopy(A); check = check) _cholesky(layout, axes, A, ::CRowMaximum; tol = 0.0, check::Bool = true) = cholesky!(cholcopy(A), CRowMaximum(); tol = tol, check = check) @@ -384,7 +384,7 @@ function _cholesky!(::SymmetricLayout{<:AbstractColumnMajor}, axes, A::AbstractM end -_inv_eye(_, ::Type{T}, axs::NTuple{2,Base.OneTo{Int}}) where T = Matrix{T}(I, map(length,axs)...) +_inv_eye(_, ::Type{T}, axs::NTuple{2,OneTo{Int}}) where T = Matrix{T}(I, map(length,axs)...) function _inv_eye(A, ::Type{T}, (rows,cols)) where T dest = zero!(similar(A, T, (cols,rows))) dest[diagind(dest)] .= one(T) @@ -414,7 +414,7 @@ macro _layoutfactorizations(Typ) LinearAlgebra.lu(A::$Typ{T}; kwds...) where T = ArrayLayouts._lu(ArrayLayouts.MemoryLayout(A), axes(A), A; kwds...) LinearAlgebra.lu!(A::$Typ, args...; kwds...) = ArrayLayouts._lu!(ArrayLayouts.MemoryLayout(A), axes(A), A, args...; kwds...) LinearAlgebra.factorize(A::$Typ) = ArrayLayouts._factorize(ArrayLayouts.MemoryLayout(A), axes(A), A) - LinearAlgebra.inv(A::$Typ) = ArrayLayouts._inv(ArrayLayouts.MemoryLayout(A), axes(A), A) + Base.inv(A::$Typ) = ArrayLayouts._inv(ArrayLayouts.MemoryLayout(A), axes(A), A) LinearAlgebra.ldiv!(L::LU{<:Any,<:$Typ}, B) = ArrayLayouts.ldiv!(L, B) LinearAlgebra.ldiv!(L::LU{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.ldiv!(L, B) end) diff --git a/src/ldiv.jl b/src/ldiv.jl index a748e41..87269bc 100644 --- a/src/ldiv.jl +++ b/src/ldiv.jl @@ -105,8 +105,8 @@ const BlasMatLdivMat{styleA, styleB, T<:BlasFloat} = MatLdivMat{styleA, styleB, const MatRdivMat{styleA, styleB, T, V} = Rdiv{styleA, styleB, <:AbstractMatrix{T}, <:AbstractMatrix{V}} const BlasMatRdivMat{styleA, styleB, T<:BlasFloat} = MatRdivMat{styleA, styleB, T, T} -materialize!(M::Ldiv{ScalarLayout}) = Base.invoke(LinearAlgebra.ldiv!, Tuple{Number,AbstractArray}, M.A, M.B) -materialize!(M::Rdiv{<:Any,ScalarLayout}) = Base.invoke(LinearAlgebra.rdiv!, Tuple{AbstractArray,Number}, M.A, M.B) +materialize!(M::Ldiv{ScalarLayout}) = invoke(LinearAlgebra.ldiv!, Tuple{Number,AbstractArray}, M.A, M.B) +materialize!(M::Rdiv{<:Any,ScalarLayout}) = invoke(LinearAlgebra.rdiv!, Tuple{AbstractArray,Number}, M.A, M.B) function materialize!(M::Ldiv{ScalarLayout,<:SymmetricLayout}) ldiv!(M.A, symmetricdata(M.B)) @@ -139,42 +139,42 @@ macro _layoutldiv(Typ) LinearAlgebra.ldiv!(A::Bidiagonal, B::$Typ) = ArrayLayouts.ldiv!(A,B) - Base.:\(A::$Typ, x::AbstractVector) = ArrayLayouts.ldiv(A,x) - Base.:\(A::$Typ, x::AbstractMatrix) = ArrayLayouts.ldiv(A,x) + (\)(A::$Typ, x::AbstractVector) = ArrayLayouts.ldiv(A,x) + (\)(A::$Typ, x::AbstractMatrix) = ArrayLayouts.ldiv(A,x) - Base.:\(x::AbstractMatrix, A::$Typ) = ArrayLayouts.ldiv(x,A) - Base.:\(x::UpperTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) - Base.:\(x::LowerTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) - Base.:\(x::Diagonal, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::AbstractMatrix, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::UpperTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::LowerTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::Diagonal, A::$Typ) = ArrayLayouts.ldiv(x,A) - Base.:\(A::Bidiagonal{<:Number}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(A,B) - Base.:\(A::Bidiagonal, B::$Typ) = ArrayLayouts.ldiv(A,B) - Base.:\(transA::Transpose{<:Number,<:Bidiagonal{<:Number}}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(transA,B) - Base.:\(transA::Transpose{<:Any,<:Bidiagonal}, B::$Typ) = ArrayLayouts.ldiv(transA,B) - Base.:\(adjA::Adjoint{<:Number,<:Bidiagonal{<:Number}}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(adjA,B) - Base.:\(adjA::Adjoint{<:Any,<:Bidiagonal}, B::$Typ) = ArrayLayouts.ldiv(adjA,B) + (\)(A::Bidiagonal{<:Number}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(A,B) + (\)(A::Bidiagonal, B::$Typ) = ArrayLayouts.ldiv(A,B) + (\)(transA::Transpose{<:Number,<:Bidiagonal{<:Number}}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(transA,B) + (\)(transA::Transpose{<:Any,<:Bidiagonal}, B::$Typ) = ArrayLayouts.ldiv(transA,B) + (\)(adjA::Adjoint{<:Number,<:Bidiagonal{<:Number}}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(adjA,B) + (\)(adjA::Adjoint{<:Any,<:Bidiagonal}, B::$Typ) = ArrayLayouts.ldiv(adjA,B) - Base.:\(x::$Typ, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::$Typ, A::$Typ) = ArrayLayouts.ldiv(x,A) - Base.:/(A::$Typ, x::AbstractVector) = ArrayLayouts.rdiv(A,x) - Base.:/(A::$Typ, x::AbstractMatrix) = ArrayLayouts.rdiv(A,x) + (/)(A::$Typ, x::AbstractVector) = ArrayLayouts.rdiv(A,x) + (/)(A::$Typ, x::AbstractMatrix) = ArrayLayouts.rdiv(A,x) - Base.:/(x::AbstractMatrix, A::$Typ) = ArrayLayouts.rdiv(x,A) - Base.:/(D::Diagonal, A::$Typ) = ArrayLayouts.rdiv(D,A) - Base.:/(A::$Typ, D::Diagonal) = ArrayLayouts.rdiv(A,D) + (/)(x::AbstractMatrix, A::$Typ) = ArrayLayouts.rdiv(x,A) + (/)(D::Diagonal, A::$Typ) = ArrayLayouts.rdiv(D,A) + (/)(A::$Typ, D::Diagonal) = ArrayLayouts.rdiv(A,D) - Base.:/(x::$Typ, A::$Typ) = ArrayLayouts.rdiv(x,A) + (/)(x::$Typ, A::$Typ) = ArrayLayouts.rdiv(x,A) end if Typ ≠ :LayoutVector ret = quote $ret - Base.:\(A::$Typ, x::LayoutVector) = ArrayLayouts.ldiv(A,x) + (\)(A::$Typ, x::LayoutVector) = ArrayLayouts.ldiv(A,x) end end if Typ ≠ :LayoutMatrix ret = quote $ret - Base.:\(A::$Typ, x::LayoutMatrix) = ArrayLayouts.ldiv(A,x) + (\)(A::$Typ, x::LayoutMatrix) = ArrayLayouts.ldiv(A,x) end end esc(ret) diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 6cfd72b..86ba62e 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -347,13 +347,13 @@ function permutelayout(layout::Union{IncreasingStrides,DecreasingStrides}, ::Val return StridedLayout() end -Base.reverse(::DenseRowMajor) = DenseColumnMajor() -Base.reverse(::RowMajor) = ColumnMajor() -Base.reverse(::DenseColumnMajor) = DenseRowMajor() -Base.reverse(::ColumnMajor) = RowMajor() -Base.reverse(::IncreasingStrides) = DecreasingStrides() -Base.reverse(::DecreasingStrides) = IncreasingStrides() -Base.reverse(::AbstractStridedLayout) = StridedLayout() +reverse(::DenseRowMajor) = DenseColumnMajor() +reverse(::RowMajor) = ColumnMajor() +reverse(::DenseColumnMajor) = DenseRowMajor() +reverse(::ColumnMajor) = RowMajor() +reverse(::IncreasingStrides) = DecreasingStrides() +reverse(::DecreasingStrides) = IncreasingStrides() +reverse(::AbstractStridedLayout) = StridedLayout() # MemoryLayout of Symmetric/Hermitian @@ -498,7 +498,7 @@ triangularlayout(::Type{Tri}, ::MemoryLayout) where Tri = Tri{UnknownLayout}() triangularlayout(::Type{Tri}, ::ML) where {Tri, ML<:AbstractColumnMajor} = Tri{ML}() triangularlayout(::Type{Tri}, ::ML) where {Tri, ML<:AbstractRowMajor} = Tri{ML}() triangularlayout(::Type{Tri}, ::ML) where {Tri, ML<:ConjLayout{<:AbstractRowMajor}} = Tri{ML}() -sublayout(layout::TriangularLayout, ::Type{<:Tuple{<:Union{Slice,Base.OneTo},<:Union{Slice,Base.OneTo}}}) = layout +sublayout(layout::TriangularLayout, ::Type{<:Tuple{<:Union{Slice,OneTo},<:Union{Slice,OneTo}}}) = layout conjlayout(::Type{<:Complex}, ::TriangularLayout{UPLO,UNIT,ML}) where {UPLO,UNIT,ML} = TriangularLayout{UPLO,UNIT,ConjLayout{ML}}() @@ -512,7 +512,7 @@ end triangulardata(A::AbstractTriangular) = parent(A) triangulardata(A::Adjoint) = Adjoint(triangulardata(parent(A))) triangulardata(A::Transpose) = Transpose(triangulardata(parent(A))) -triangulardata(A::SubArray{<:Any,2,<:Any,<:Tuple{<:Union{Slice,Base.OneTo},<:Union{Slice,Base.OneTo}}}) = +triangulardata(A::SubArray{<:Any,2,<:Any,<:Tuple{<:Union{Slice,OneTo},<:Union{Slice,OneTo}}}) = view(triangulardata(parent(A)), parentindices(A)...) ### diff --git a/src/mul.jl b/src/mul.jl index dfbfc3b..23672ca 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -121,23 +121,23 @@ broadcastable(M::Mul) = M macro veclayoutmul(Typ) ret = quote - Base.:*(A::AbstractMatrix, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::Adjoint{<:Any,<:AbstractMatrix{T}}, B::$Typ{S}) where {T,S} = ArrayLayouts.mul(A,B) - Base.:*(A::Transpose{<:Any,<:AbstractMatrix{T}}, B::$Typ{S}) where {T,S} = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AdjointAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AdjointAbsVec{<:Number}, B::$Typ{<:Number}) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec{T}, B::$Typ{T}) where T<:Real = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec{T,<:$Typ{T}}, B::$Typ{T}) where T<:Real = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec{T,<:$Typ{T}}, B::AbstractVector{T}) where T<:Real = ArrayLayouts.mul(A,B) - - Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::AbstractMatrix, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::Adjoint{<:Any,<:AbstractMatrix{T}}, B::$Typ{S}) where {T,S} = ArrayLayouts.mul(A,B) + (*)(A::Transpose{<:Any,<:AbstractMatrix{T}}, B::$Typ{S}) where {T,S} = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointAbsVec{<:Number}, B::$Typ{<:Number}) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec{T}, B::$Typ{T}) where T<:Real = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec{T,<:$Typ{T}}, B::$Typ{T}) where T<:Real = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec{T,<:$Typ{T}}, B::AbstractVector{T}) where T<:Real = ArrayLayouts.mul(A,B) + + (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) end for Struc in (:AbstractTriangular, :Diagonal) ret = quote $ret - Base.:*(A::LinearAlgebra.$Struc, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.$Struc, B::$Typ) = ArrayLayouts.mul(A,B) end end for Mod in (:AdjointAbsVec, :TransposeAbsVec) @@ -147,14 +147,14 @@ macro veclayoutmul(Typ) LinearAlgebra.mul!(dest::AbstractVector, A::$Mod{<:Any,<:$Typ}, b::AbstractVector) = ArrayLayouts.mul!(dest,A,b) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::ArrayLayouts.LayoutMatrix) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::AbstractMatrix) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::AbstractVector) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Number,<:$Typ}, B::AbstractVector{<:Number}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Number,<:$Typ}, B::$Typ{<:Number}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::Diagonal) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::LinearAlgebra.AbstractTriangular) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::ArrayLayouts.LayoutMatrix) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::AbstractMatrix) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::AbstractVector) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Number,<:$Typ}, B::AbstractVector{<:Number}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Number,<:$Typ}, B::$Typ{<:Number}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::Diagonal) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::LinearAlgebra.AbstractTriangular) = ArrayLayouts.mul(A,B) end end @@ -183,35 +183,35 @@ macro layoutmul(Typ) LinearAlgebra.mul!(dest::AbstractMatrix, A::$Typ, B::$Typ, α::Number, β::Number) = ArrayLayouts.mul!(dest,A,B,α,β) - Base.:*(A::$Typ, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::AbstractMatrix) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::AbstractVector) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::ArrayLayouts.LayoutVector) = ArrayLayouts.mul(A,B) - Base.:*(A::AbstractMatrix, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AdjointAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AdjointAbsVec{<:Any,<:Zeros{<:Any,1}}, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec{<:Any,<:Zeros{<:Any,1}}, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.Transpose{T,<:$Typ}, B::Zeros{T,1}) where T<:Real = ArrayLayouts.mul(A,B) - - Base.:*(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::AbstractMatrix) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::AbstractVector) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::ArrayLayouts.LayoutVector) = ArrayLayouts.mul(A,B) + (*)(A::AbstractMatrix, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointAbsVec{<:Any,<:Zeros{<:Any,1}}, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec{<:Any,<:Zeros{<:Any,1}}, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.Transpose{T,<:$Typ}, B::Zeros{T,1}) where T<:Real = ArrayLayouts.mul(A,B) + + (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) # disambiguation - Base.:*(A::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}, B::$Typ) = + (*)(A::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}) = + (*)(A::$Typ, B::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}) = ArrayLayouts.mul(A,B) if isdefined(LinearAlgebra, :AdjointQ) - Base.:*(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) end end for Struc in (:AbstractTriangular, :Diagonal) ret = quote $ret - Base.:*(A::LinearAlgebra.$Struc, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::LinearAlgebra.$Struc) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.$Struc, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::LinearAlgebra.$Struc) = ArrayLayouts.mul(A,B) end end for Mod in (:Adjoint, :Transpose, :Symmetric, :Hermitian) @@ -237,23 +237,23 @@ macro layoutmul(Typ) LinearAlgebra.mul!(dest::AbstractMatrix, A::$Mod{<:Any,<:$Typ}, B::$Mod{<:Any,<:$Typ}, α::Number, β::Number) = ArrayLayouts.mul!(dest,A,B,α,β) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::AbstractMatrix) = ArrayLayouts.mul(A,B) - Base.:*(A::AbstractMatrix, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AdjointAbsVec, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.TransposeAbsVec, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::AbstractVector) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::ArrayLayouts.LayoutVector) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::Zeros{<:Any,1}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::AbstractMatrix) = ArrayLayouts.mul(A,B) + (*)(A::AbstractMatrix, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointAbsVec, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.TransposeAbsVec, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::AbstractVector) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::ArrayLayouts.LayoutVector) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::Zeros{<:Any,1}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.mul(A,B) - Base.:*(A::$Typ, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::Diagonal) = ArrayLayouts.mul(A,B) - Base.:*(A::Diagonal, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::Diagonal) = ArrayLayouts.mul(A,B) + (*)(A::Diagonal, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::LinearAlgebra.AbstractTriangular, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) - Base.:*(A::$Mod{<:Any,<:$Typ}, B::LinearAlgebra.AbstractTriangular) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AbstractTriangular, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B) + (*)(A::$Mod{<:Any,<:$Typ}, B::LinearAlgebra.AbstractTriangular) = ArrayLayouts.mul(A,B) end end @@ -282,7 +282,7 @@ struct Dot{StyleA,StyleB,ATyp,BTyp} end @inline Dot(A::ATyp,B::BTyp) where {ATyp,BTyp} = Dot{typeof(MemoryLayout(ATyp)), typeof(MemoryLayout(BTyp)), ATyp, BTyp}(A, B) -@inline copy(d::Dot) = Base.invoke(LinearAlgebra.dot, Tuple{AbstractArray,AbstractArray}, d.A, d.B) +@inline copy(d::Dot) = invoke(LinearAlgebra.dot, Tuple{AbstractArray,AbstractArray}, d.A, d.B) @inline materialize(d::Dot) = copy(instantiate(d)) @inline Dot(M::Mul{<:DualLayout,<:Any,<:AbstractMatrix,<:AbstractVector}) = Dot(M.A', M.B) @inline mulreduce(M::Mul{<:DualLayout,<:Any,<:AbstractMatrix,<:AbstractVector}) = Dot(M) @@ -312,7 +312,7 @@ LinearAlgebra.dot(x::AbstractVector, A::Symmetric{<:Real,<:LayoutMatrix}, y::Abs # allow overloading for infinite or lazy case -@inline _power_by_squaring(_, _, A, p) = Base.invoke(Base.power_by_squaring, Tuple{AbstractMatrix,Integer}, A, p) +@inline _power_by_squaring(_, _, A, p) = invoke(Base.power_by_squaring, Tuple{AbstractMatrix,Integer}, A, p) # TODO: Remove unnecessary _apply _apply(_, _, op, Λ::UniformScaling, A::AbstractMatrix) = op(Diagonal(Fill(Λ.λ,(axes(A,1),))), A) _apply(_, _, op, A::AbstractMatrix, Λ::UniformScaling) = op(A, Diagonal(Fill(Λ.λ,(axes(A,1),)))) diff --git a/src/muladd.jl b/src/muladd.jl index c1a9621..a171a40 100644 --- a/src/muladd.jl +++ b/src/muladd.jl @@ -221,7 +221,7 @@ function _default_blasmul!(::IndexCartesian, α, A::AbstractMatrix, B::AbstractV end default_blasmul!(α, A::AbstractMatrix, B::AbstractVector, β, C::AbstractVector) = - _default_blasmul!(Base.IndexStyle(typeof(A)), α, A, B, β, C) + _default_blasmul!(IndexStyle(typeof(A)), α, A, B, β, C) function materialize!(M::MatMulMatAdd) α, A, B, β, C = M.α, M.A, M.B, M.β, M.C From 5c0788dcec62c0b9fe9e17d11c0e3affdab20e23 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sat, 5 Nov 2022 17:39:45 +0100 Subject: [PATCH 11/16] rm AdjointQ stuff --- src/ArrayLayouts.jl | 2 -- src/factorizations.jl | 32 ++++++++++++++++---------------- src/lmul.jl | 4 ++-- src/memorylayout.jl | 3 --- src/mul.jl | 10 +--------- 5 files changed, 19 insertions(+), 32 deletions(-) diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index d8117c5..879d9e3 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -21,8 +21,6 @@ using LinearAlgebra: AbstractTriangular, AbstractQ, QRCompactWYQ, QRPackedQ, che AdjOrTrans, HermOrSym, RealHermSymComplexHerm, AdjointAbsVec, TransposeAbsVec, checknonsingular, _apply_ipiv_rows!, ipiv2perm, chkfullrank -AdjointQtype{T} = isdefined(LinearAlgebra, :AdjointQ) ? LinearAlgebra.AdjointQ{T} : Adjoint{T,<:AbstractQ} - using LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex using FillArrays: AbstractFill, getindex_value, axes_print_matrix_row, _copy_oftype diff --git a/src/factorizations.jl b/src/factorizations.jl index 26608fb..3280226 100644 --- a/src/factorizations.jl +++ b/src/factorizations.jl @@ -196,14 +196,14 @@ end ### QcB -materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = - (A = parent(M.A); LAPACK.ormqr!('L','T',A.factors,A.τ,M.B)) -materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = - (A = parent(M.A); LAPACK.ormqr!('L','C',A.factors,A.τ,M.B)) -materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = - (A = parent(M.A); LAPACK.gemqrt!('L','T',A.factors,A.T,M.B)) -materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AdjointQtype{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = - (A = parent(M.A); LAPACK.gemqrt!('L','C',A.factors,A.T,M.B)) +materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = + (A = M.A.parent; LAPACK.ormqr!('L','T',A.factors,A.τ,M.B)) +materialize!(M::Lmul{<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = + (A = M.A.parent; LAPACK.ormqr!('L','C',A.factors,A.τ,M.B)) +materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasFloat = + (A = M.A.parent; LAPACK.gemqrt!('L','T',A.factors,A.T,M.B)) +materialize!(M::Lmul{<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractStridedLayout,<:AbstractMatrix{T},<:AbstractVecOrMat{T}}) where T<:BlasComplex = + (A = M.A.parent; LAPACK.gemqrt!('L','C',A.factors,A.T,M.B)) function materialize!(M::Lmul{<:AdjQRPackedQLayout}) adjA,B = M.A, M.B require_one_based_indexing(B) @@ -264,14 +264,14 @@ function materialize!(M::Rmul{<:Any,<:QRPackedQLayout}) end ### AQc -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasReal = - (B = parent(M.B); LAPACK.ormqr!('R','T',B.factors,B.τ,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasComplex = - (B = parent(M.B); LAPACK.ormqr!('R','C',B.factors,B.τ,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasReal = - (B = parent(M.B); LAPACK.gemqrt!('R','T',B.factors,B.T,M.A)) -materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AdjointQtype{T}}) where T<:BlasComplex = - (B = parent(M.B); LAPACK.gemqrt!('R','C',B.factors,B.T,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasReal = + (B = M.B.parent; LAPACK.ormqr!('R','T',B.factors,B.τ,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRPackedQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasComplex = + (B = M.B.parent; LAPACK.ormqr!('R','C',B.factors,B.τ,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasReal = + (B = M.B.parent; LAPACK.gemqrt!('R','T',B.factors,B.T,M.A)) +materialize!(M::Rmul{<:AbstractStridedLayout,<:AdjQRCompactWYQLayout{<:AbstractStridedLayout,<:AbstractStridedLayout},<:AbstractVecOrMat{T},<:AbstractMatrix{T}}) where T<:BlasComplex = + (B = M.B.parent; LAPACK.gemqrt!('R','C',B.factors,B.T,M.A)) function materialize!(M::Rmul{<:Any,<:AdjQRPackedQLayout}) A,adjQ = M.A,M.B Q = parent(adjQ) diff --git a/src/lmul.jl b/src/lmul.jl index 253062d..6ce9818 100644 --- a/src/lmul.jl +++ b/src/lmul.jl @@ -70,8 +70,8 @@ copyto!(dest::AbstractArray, M::Rmul) = _rmul_copyto!(dest, M) materialize!(M::Lmul) = LinearAlgebra.lmul!(M.A,M.B) materialize!(M::Rmul) = LinearAlgebra.rmul!(M.A,M.B) -materialize!(M::Lmul{ScalarLayout}) = Base.invoke(LinearAlgebra.lmul!, Tuple{Number,AbstractArray}, M.A, M.B) -materialize!(M::Rmul{<:Any,ScalarLayout}) = Base.invoke(LinearAlgebra.rmul!, Tuple{AbstractArray,Number}, M.A, M.B) +materialize!(M::Lmul{ScalarLayout}) = invoke(LinearAlgebra.lmul!, Tuple{Number,AbstractArray}, M.A, M.B) +materialize!(M::Rmul{<:Any,ScalarLayout}) = invoke(LinearAlgebra.rmul!, Tuple{AbstractArray,Number}, M.A, M.B) function materialize!(M::Lmul{ScalarLayout,<:SymmetricLayout}) lmul!(M.A, symmetricdata(M.B)) diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 86ba62e..c7fbe33 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -281,9 +281,6 @@ struct DualLayout{ML<:MemoryLayout} <: MemoryLayout end MemoryLayout(::Type{Transpose{T,P}}) where {T,P} = transposelayout(MemoryLayout(P)) MemoryLayout(::Type{Adjoint{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) -if isdefined(LinearAlgebra, :AdjointQ) - MemoryLayout(::Type{LinearAlgebra.AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) -end MemoryLayout(::Type{AdjointAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(adjointlayout(T,MemoryLayout(P)))}() MemoryLayout(::Type{TransposeAbsVec{T,P}}) where {T,P<:AbstractVector} = DualLayout{typeof(transposelayout(MemoryLayout(P)))}() diff --git a/src/mul.jl b/src/mul.jl index 23672ca..4ac1f88 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -132,6 +132,7 @@ macro veclayoutmul(Typ) (*)(A::LinearAlgebra.TransposeAbsVec{T,<:$Typ{T}}, B::AbstractVector{T}) where T<:Real = ArrayLayouts.mul(A,B) (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::LinearAlgebra.LQPackedQ) = ArrayLayouts.mul(A,B) end for Struc in (:AbstractTriangular, :Diagonal) ret = quote @@ -196,15 +197,6 @@ macro layoutmul(Typ) (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) (*)(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) - # disambiguation - (*)(A::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}, B::$Typ) = - ArrayLayouts.mul(A,B) - (*)(A::$Typ, B::Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ}) = - ArrayLayouts.mul(A,B) - if isdefined(LinearAlgebra, :AdjointQ) - (*)(A::LinearAlgebra.AdjointQ, B::$Typ) = ArrayLayouts.mul(A,B) - (*)(A::$Typ, B::LinearAlgebra.AdjointQ) = ArrayLayouts.mul(A,B) - end end for Struc in (:AbstractTriangular, :Diagonal) ret = quote From bc3fb7d5c97e41c4dff0521759e1a1a3d2ccde90 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sun, 6 Nov 2022 08:36:54 +0100 Subject: [PATCH 12/16] bump patch version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 3c0fb86..1e54175 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ArrayLayouts" uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" authors = ["Sheehan Olver "] -version = "0.8.12" +version = "0.8.13" [deps] FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" From 15130826c8c9532458aeb3bc311462fef5a5dc93 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sun, 6 Nov 2022 10:58:08 +0100 Subject: [PATCH 13/16] add tests and fix ambiguities --- src/ldiv.jl | 2 ++ test/test_layoutarray.jl | 13 +++++++++---- 2 files changed, 11 insertions(+), 4 deletions(-) diff --git a/src/ldiv.jl b/src/ldiv.jl index 87269bc..73ada2d 100644 --- a/src/ldiv.jl +++ b/src/ldiv.jl @@ -144,7 +144,9 @@ macro _layoutldiv(Typ) (\)(x::AbstractMatrix, A::$Typ) = ArrayLayouts.ldiv(x,A) (\)(x::UpperTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::UnitUpperTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) (\)(x::LowerTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) + (\)(x::UnitLowerTriangular, A::$Typ) = ArrayLayouts.ldiv(x,A) (\)(x::Diagonal, A::$Typ) = ArrayLayouts.ldiv(x,A) (\)(A::Bidiagonal{<:Number}, B::$Typ{<:Number}) = ArrayLayouts.ldiv(A,B) diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index dd78487..816ef0f 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -66,8 +66,8 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() end b = randn(5) for Tri in (UpperTriangular, UnitUpperTriangular, LowerTriangular, UnitLowerTriangular) - @test ldiv!(Tri(A), copy(b)) ≈ ldiv!(Tri(A.A), copy(b)) - @test lmul!(Tri(A), copy(b)) ≈ lmul!(Tri(A.A), copy(b)) + @test ldiv!(Tri(A), copy(b)) ≈ ldiv!(Tri(A.A), copy(b)) ≈ Tri(A.A) \ MyVector(b) + @test lmul!(Tri(A), copy(b)) ≈ lmul!(Tri(A.A), copy(b)) ≈ Tri(A.A) * MyVector(b) end @test copyto!(MyMatrix(Array{Float64}(undef,5,5)), A) == A @@ -94,7 +94,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() @test_throws ErrorException qr!(A) @test lu!(copy(A)).factors ≈ lu(A.A).factors b = randn(5) - @test all(A \ b .≡ A.A \ b) + @test all(A \ b .≡ A.A \ b .≡ A.A \ MyVector(b)) @test all(lu(A).L .≡ lu(A.A).L) @test all(lu(A).U .≡ lu(A.A).U) @test lu(A).p == lu(A.A).p @@ -105,7 +105,10 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() S = Symmetric(MyMatrix(reshape(inv.(1:25),5,5) + 10I)) @test cholesky(S).U ≈ @inferred(cholesky!(deepcopy(S))).U - @test cholesky(S,CRowMaximum()).U ≈ cholesky(Matrix(S),CRowMaximum()).U + @test cholesky(S, CRowMaximum()).U ≈ cholesky(Matrix(S), CRowMaximum()).U + @test cholesky(S) \ b ≈ cholesky(Matrix(S)) \ b ≈ cholesky(Matrix(S)) \ MyVector(b) + @test cholesky(S, CRowMaximum()) \ b ≈ cholesky(Matrix(S), CRowMaximum()) \ b + @test cholesky(S, CRowMaximum()) \ b ≈ cholesky(Matrix(S), CRowMaximum()) \ MyVector(b) S = Symmetric(MyMatrix(reshape(inv.(1:25),5,5) + 10I),:L) @test cholesky(S).U ≈ @inferred(cholesky!(deepcopy(S))).U @@ -164,6 +167,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() t̃ = copy(t) T̃ = copy(T) B = Bidiagonal(randn(5),randn(4),:U) + D = Diagonal(randn(5)) @test ldiv!(A, copy(x)) ≈ A\x @test A\t ≈ A\t̃ # QR is not general enough @@ -173,6 +177,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() @test_broken A/T ≈ A/T̃ @test_broken ldiv!(A, T) ≈ A\T @test B\A ≈ B\Matrix(A) + @test D \ A ≈ D \ Matrix(A) @test transpose(B)\A ≈ transpose(B)\Matrix(A) ≈ Transpose(B)\A ≈ Adjoint(B)\A @test B'\A ≈ B'\Matrix(A) @test A\A ≈ I From c3cfd3f631bb12d4aed61ee2fabc215013443c0d Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sun, 6 Nov 2022 11:08:27 +0100 Subject: [PATCH 14/16] run CI on nightly --- .github/workflows/ci.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 1cc4c47..09d6741 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -14,6 +14,7 @@ jobs: version: - '1.6' - '1' + - 'nightly' os: - ubuntu-latest - macOS-latest From c5d1d3163ae1bb0d6596602d20f11c9bd2d4b55d Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sun, 6 Nov 2022 11:51:16 +0100 Subject: [PATCH 15/16] add factorization tests --- test/test_layoutarray.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index e4eb62c..fc20cbe 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -93,7 +93,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() @test_throws ErrorException qr!(A) @test lu!(copy(A)).factors ≈ lu(A.A).factors b = randn(5) - @test all(A \ b .≡ A.A \ b .≡ A.A \ MyVector(b)) + @test all(A \ b .≡ A.A \ b .≡ A.A \ MyVector(b) .≡ ldiv!(lu(A.A)), copy(MyVector(b))) @test all(lu(A).L .≡ lu(A.A).L) @test all(lu(A).U .≡ lu(A.A).U) @test lu(A).p == lu(A.A).p @@ -107,7 +107,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() @test cholesky(S, CRowMaximum()).U ≈ cholesky(Matrix(S), CRowMaximum()).U @test cholesky(S) \ b ≈ cholesky(Matrix(S)) \ b ≈ cholesky(Matrix(S)) \ MyVector(b) @test cholesky(S, CRowMaximum()) \ b ≈ cholesky(Matrix(S), CRowMaximum()) \ b - @test cholesky(S, CRowMaximum()) \ b ≈ cholesky(Matrix(S), CRowMaximum()) \ MyVector(b) + @test cholesky(S, CRowMaximum()) \ b ≈ ldiv!(cholesky(Matrix(S), CRowMaximum()), copy(MyVector(b))) S = Symmetric(MyMatrix(reshape(inv.(1:25),5,5) + 10I),:L) @test cholesky(S).U ≈ @inferred(cholesky!(deepcopy(S))).U From 56ee5102bd1f1cb80d10803d56050c0a231dd114 Mon Sep 17 00:00:00 2001 From: Daniel Karrasch Date: Sun, 6 Nov 2022 12:00:49 +0100 Subject: [PATCH 16/16] fix brackets --- test/test_layoutarray.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index fc20cbe..9360bb7 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -93,7 +93,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor() @test_throws ErrorException qr!(A) @test lu!(copy(A)).factors ≈ lu(A.A).factors b = randn(5) - @test all(A \ b .≡ A.A \ b .≡ A.A \ MyVector(b) .≡ ldiv!(lu(A.A)), copy(MyVector(b))) + @test all(A \ b .≡ A.A \ b .≡ A.A \ MyVector(b) .≡ ldiv!(lu(A.A), copy(MyVector(b)))) @test all(lu(A).L .≡ lu(A.A).L) @test all(lu(A).U .≡ lu(A.A).U) @test lu(A).p == lu(A.A).p