diff --git a/Project.toml b/Project.toml index b77f3b8..2b80ef4 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ArrayLayouts" uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" -version = "1.12.2" +version = "1.12.3" authors = ["Sheehan Olver "] [deps] diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl index 6aade6a..738ed1b 100644 --- a/src/ArrayLayouts.jl +++ b/src/ArrayLayouts.jl @@ -3,7 +3,7 @@ using Base: _typed_hcat using Base.Broadcast, LinearAlgebra, FillArrays using LinearAlgebra.BLAS -using Base: AbstractCartesianIndex, OneTo, oneto, RangeIndex, ReinterpretArray, ReshapedArray, +using Base: AbstractCartesianIndex, LazyString, OneTo, oneto, RangeIndex, ReinterpretArray, ReshapedArray, Slice, tuple_type_tail, unalias, @propagate_inbounds @@ -21,7 +21,7 @@ using LinearAlgebra: AbstractQ, AbstractTriangular, AdjOrTrans, AdjointAbsVec, H using LinearAlgebra.BLAS: BlasComplex, BlasFloat, BlasReal -const AdjointQtype{T} = isdefined(LinearAlgebra, :AdjointQ) ? LinearAlgebra.AdjointQ{T} : Adjoint{T,<:AbstractQ} +const AdjointQtype{T} = LinearAlgebra.AdjointQ{T} using FillArrays: AbstractFill, axes_print_matrix_row, getindex_value using StaticArrays @@ -37,19 +37,9 @@ export materialize, materialize!, MulAdd, muladd!, Ldiv, Rdiv, Lmul, Rmul, Dot, colsupport, rowsupport, layout_getindex, AbstractQLayout, LayoutArray, LayoutMatrix, LayoutVector, RangeCumsum -if VERSION < v"1.7-" - const ColumnNorm = Val{true} - const RowMaximum = Val{true} - const NoPivot = Val{false} -end - -if VERSION < v"1.8-" - const CRowMaximum = Val{true} - const CNoPivot = Val{false} -else - const CRowMaximum = RowMaximum - const CNoPivot = NoPivot -end +# retained for backwards compatibility with downstream packages +const CRowMaximum = RowMaximum +const CNoPivot = NoPivot if VERSION ≥ v"1.11.0-DEV.21" using LinearAlgebra: UpperOrLowerTriangular @@ -60,12 +50,6 @@ else LinearAlgebra.UnitLowerTriangular{T,S}} end -@static if VERSION ≥ v"1.8.0" - import Base: LazyString -else - const LazyString = string -end - # Originally defined in FillArrays _copy_oftype(A::AbstractArray, ::Type{S}) where {S} = eltype(A) == S ? copy(A) : AbstractArray{S}(A) _copy_oftype(A::AbstractRange, ::Type{S}) where {S} = eltype(A) == S ? copy(A) : map(S, A) @@ -200,6 +184,11 @@ for Typ in (:LayoutArray, :(Transpose{<:Any,<:LayoutMatrix}), :(Adjoint{<:Any,<: end end +# disambiguation with the LinearAlgebra methods for transposes of real/complex +# matrices, which were narrowed to `Transpose{<:Union{Real,Complex}}` in Julia v1.13 +LinearAlgebra.lmul!(α::Union{Real,Complex}, A::Transpose{<:Union{Real,Complex},<:LayoutMatrix}) = lmul!(α, A) +LinearAlgebra.rmul!(A::Transpose{<:Union{Real,Complex},<:LayoutMatrix}, α::Union{Real,Complex}) = rmul!(A, α) + getindex(A::LayoutVector, kr::AbstractVector) = layout_getindex(A, kr) getindex(A::LayoutVector, kr::Colon) = layout_getindex(A, kr) getindex(A::AdjOrTrans{<:Any,<:LayoutVector}, kr::Integer, jr::Colon) = layout_getindex(A, kr, jr) @@ -287,11 +276,9 @@ copyto!(dest::AbstractMatrix, src::AdjOrTrans{<:Any,<:LayoutArray}) = copyto!_la copyto!(dest::SubArray{<:Any,2,<:LayoutArray}, src::AdjOrTrans{<:Any,<:LayoutArray}) = copyto!_layout(dest, src) copyto!(dest::SubArray{<:Any,2,<:LayoutMatrix}, src::SubArray{<:Any,2,<:AdjOrTrans{<:Any,<:LayoutArray}}) = copyto!_layout(dest, src) copyto!(dest::AbstractMatrix, src::SubArray{<:Any,2,<:AdjOrTrans{<:Any,<:LayoutArray}}) = copyto!_layout(dest, src) -if isdefined(LinearAlgebra, :copymutable_oftype) - LinearAlgebra.copymutable_oftype(A::Union{LayoutArray,Symmetric{<:Any,<:LayoutMatrix},Hermitian{<:Any,<:LayoutMatrix}, - UpperOrLowerTriangular{<:Any,<:LayoutMatrix}, - AdjOrTrans{<:Any,<:LayoutMatrix}}, ::Type{T}) where T = copymutable_oftype_layout(MemoryLayout(A), A, T) -end +LinearAlgebra.copymutable_oftype(A::Union{LayoutArray,Symmetric{<:Any,<:LayoutMatrix},Hermitian{<:Any,<:LayoutMatrix}, + UpperOrLowerTriangular{<:Any,<:LayoutMatrix}, + AdjOrTrans{<:Any,<:LayoutMatrix}}, ::Type{T}) where T = copymutable_oftype_layout(MemoryLayout(A), A, T) copymutable_oftype_layout(_, A, ::Type{S}) where S = copyto!(similar(A, S), A) # avoid bad copy in Base diff --git a/src/factorizations.jl b/src/factorizations.jl index 2cfc3cd..6e09a13 100644 --- a/src/factorizations.jl +++ b/src/factorizations.jl @@ -325,26 +325,14 @@ factorize_layout(layout, axes, A) = qr(A) # Default to QR factorize_layout(::AbstractStridedLayout, axes, A) = lu(A) -if VERSION < v"1.8-" - function lu!_layout(::AbstractColumnMajor, axes, A::AbstractMatrix{T}, pivot::Union{NoPivot, RowMaximum} = RowMaximum(); - check::Bool = true) where T<:BlasFloat - if pivot === NoPivot() - return generic_lufact!(A, pivot; check = check) - end - lpt = LAPACK.getrf!(A) - check && checknonsingular(lpt[3]) - return LU{T,typeof(A)}(lpt[1], lpt[2], lpt[3]) - end -else - function lu!_layout(::AbstractColumnMajor, axes, A::AbstractMatrix{T}, pivot::Union{NoPivot, RowMaximum} = RowMaximum(); - check::Bool = true) where T<:BlasFloat - if pivot === NoPivot() - return generic_lufact!(A, pivot; check = check) - end - lpt = LAPACK.getrf!(A) - check && checknonsingular(lpt[3]) - return LU{T,typeof(A),typeof(lpt[2])}(lpt[1], lpt[2], lpt[3]) +function lu!_layout(::AbstractColumnMajor, axes, A::AbstractMatrix{T}, pivot::Union{NoPivot, RowMaximum} = RowMaximum(); + check::Bool = true) where T<:BlasFloat + if pivot === NoPivot() + return generic_lufact!(A, pivot; check = check) end + lpt = LAPACK.getrf!(A) + check && checknonsingular(lpt[3]) + return LU{T,typeof(A),typeof(lpt[2])}(lpt[1], lpt[2], lpt[3]) end # for some reason only defined for StridedMatrix in LinearAlgebra @@ -390,11 +378,7 @@ end function cholesky!_layout(::SymmetricLayout{<:AbstractColumnMajor}, axes, A::AbstractMatrix{<:BlasReal}, ::CRowMaximum; tol = 0.0, check::Bool = true) AA, piv, rank, info = LAPACK.pstrf!(A.uplo, A.data, tol) - if VERSION < v"1.8" - C = CholeskyPivoted{eltype(AA),typeof(AA)}(AA, A.uplo, piv, rank, tol, info) - else - C = CholeskyPivoted{eltype(AA),typeof(AA),typeof(piv)}(AA, A.uplo, piv, rank, tol, info) - end + C = CholeskyPivoted{eltype(AA),typeof(AA),typeof(piv)}(AA, A.uplo, piv, rank, tol, info) check && chkfullrank(C) return C end diff --git a/src/ldiv.jl b/src/ldiv.jl index 23d076c..689463f 100644 --- a/src/ldiv.jl +++ b/src/ldiv.jl @@ -83,13 +83,7 @@ __ldiv!(_, F, B) = LinearAlgebra.ldiv!(F, B) @inline _ldiv!(dest, A, B; kwds...) = ldiv!(dest, factorize(A), B; kwds...) @inline _ldiv!(dest, A::Factorization, B; kwds...) = LinearAlgebra.ldiv!(dest, A, B; kwds...) -if VERSION ≥ v"1.10-" - using LinearAlgebra: TransposeFactorization, AdjointFactorization -else - const TransposeFactorization = Transpose - const AdjointFactorization = Adjoint - -end +using LinearAlgebra: TransposeFactorization, AdjointFactorization @inline _ldiv!(dest, A::TransposeFactorization{<:Any,<:Factorization}, B; kwds...) = LinearAlgebra.ldiv!(dest, A, B; kwds...) @inline _ldiv!(dest, A::AdjointFactorization{<:Any,<:Factorization}, B; kwds...) = LinearAlgebra.ldiv!(dest, A, B; kwds...) @@ -164,9 +158,6 @@ macro _layoutldiv(Typ) (\)(x::AbstractMatrix, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) (\)(x::LinearAlgebra.HermOrSym, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) - if VERSION < v"1.9-" # disambiguation - \(x::LinearAlgebra.HermOrSym{<:Any,<:StridedMatrix}, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) - end (\)(x::UpperTriangular, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) (\)(x::UnitUpperTriangular, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) (\)(x::LowerTriangular, A::$Typ; kwds...) = ArrayLayouts.ldiv(x,A; kwds...) diff --git a/src/memorylayout.jl b/src/memorylayout.jl index 2570a11..a887a48 100644 --- a/src/memorylayout.jl +++ b/src/memorylayout.jl @@ -283,9 +283,7 @@ MemoryLayout(::Type{Transpose{T,P}}) where {T,P} = transposelayout(MemoryLayout( MemoryLayout(::Type{Adjoint{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) 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)))}() -if isdefined(LinearAlgebra, :AdjointQ) - MemoryLayout(::Type{LinearAlgebra.AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) -end +MemoryLayout(::Type{LinearAlgebra.AdjointQ{T,P}}) where {T,P} = adjointlayout(T, MemoryLayout(P)) transposelayout(_) = UnknownLayout() transposelayout(::StridedLayout) = StridedLayout() diff --git a/src/mul.jl b/src/mul.jl index 377c397..0cb89cf 100644 --- a/src/mul.jl +++ b/src/mul.jl @@ -161,16 +161,14 @@ macro veclayoutmul(Typ) (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) (*)(A::$Typ, B::LinearAlgebra.LQPackedQ) = ArrayLayouts.mul(A,B) end - if isdefined(LinearAlgebra, :AdjointQ) - ret = quote - $ret + ret = quote + $ret - const FlexibleLeftQs = Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ} - # disambiguation for flexible left-mul Qs - (*)(A::FlexibleLeftQs, B::$Typ) = ArrayLayouts.mul(A,B) - # flexible right-mul/adjoint left-mul Qs - (*)(A::LinearAlgebra.AdjointQ{<:Any,<:LinearAlgebra.LQPackedQ}, B::$Typ) = ArrayLayouts.mul(A,B) - end + const FlexibleLeftQs = Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ} + # disambiguation for flexible left-mul Qs + (*)(A::FlexibleLeftQs, B::$Typ) = ArrayLayouts.mul(A,B) + # flexible right-mul/adjoint left-mul Qs + (*)(A::LinearAlgebra.AdjointQ{<:Any,<:LinearAlgebra.LQPackedQ}, B::$Typ) = ArrayLayouts.mul(A,B) end for Struc in (:AbstractTriangular, :Diagonal) ret = quote @@ -234,18 +232,16 @@ macro layoutmul(Typ) (*)(A::LinearAlgebra.AbstractQ, B::$Typ) = ArrayLayouts.mul(A,B) (*)(A::$Typ, B::LinearAlgebra.AbstractQ) = ArrayLayouts.mul(A,B) end - if isdefined(LinearAlgebra, :AdjointQ) - ret = quote - $ret + ret = quote + $ret - const FlexibleLeftQs = Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ} - # disambiguation for flexible left-mul/adjoint right-mul Qs - (*)(A::FlexibleLeftQs, B::$Typ) = ArrayLayouts.mul(A,B) - (*)(A::$Typ, B::LinearAlgebra.AdjointQ{<:Any,<:FlexibleLeftQs}) = ArrayLayouts.mul(A,B) - # disambiguation for flexible right-mul/adjoint left-mul Qs - (*)(A::$Typ, B::LinearAlgebra.LQPackedQ) = ArrayLayouts.mul(A,B) - (*)(A::LinearAlgebra.AdjointQ{<:Any,<:LinearAlgebra.LQPackedQ}, B::$Typ) = ArrayLayouts.mul(A,B) - end + const FlexibleLeftQs = Union{LinearAlgebra.HessenbergQ, LinearAlgebra.QRCompactWYQ, LinearAlgebra.QRPackedQ} + # disambiguation for flexible left-mul/adjoint right-mul Qs + (*)(A::FlexibleLeftQs, B::$Typ) = ArrayLayouts.mul(A,B) + (*)(A::$Typ, B::LinearAlgebra.AdjointQ{<:Any,<:FlexibleLeftQs}) = ArrayLayouts.mul(A,B) + # disambiguation for flexible right-mul/adjoint left-mul Qs + (*)(A::$Typ, B::LinearAlgebra.LQPackedQ) = ArrayLayouts.mul(A,B) + (*)(A::LinearAlgebra.AdjointQ{<:Any,<:LinearAlgebra.LQPackedQ}, B::$Typ) = ArrayLayouts.mul(A,B) end for Struc in (:AbstractTriangular, :Diagonal, :Bidiagonal, :SymTridiagonal, :Tridiagonal) # starting from Julia v1.10, the last four could be put into a single Union to diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl index cd8f4f3..09e7543 100644 --- a/test/test_layoutarray.jl +++ b/test/test_layoutarray.jl @@ -104,10 +104,8 @@ BroadcastStyle(::Type{<:MyVector2{T}}) where {T} = MyBroadcastStyle{1}() for Tri in (UpperTriangular, UnitUpperTriangular, LowerTriangular, UnitLowerTriangular) @test ldiv!(Tri(A), copy(b)) ≈ ldiv!(Tri(A.A), copy(b)) ≈ Tri(A.A) \ MyVector(b) @test ldiv!(Tri(A), copy(B)) ≈ ldiv!(Tri(A.A), copy(B)) ≈ Tri(A.A) \ MyMatrix(B) - if VERSION ≥ v"1.9" - @test rdiv!(copy(b)', Tri(A)) ≈ rdiv!(copy(b)', Tri(A.A)) ≈ MyVector(b)' / Tri(A.A) - @test rdiv!(copy(B), Tri(A)) ≈ rdiv!(copy(B), Tri(A.A)) ≈ B / Tri(A.A) - end + @test rdiv!(copy(b)', Tri(A)) ≈ rdiv!(copy(b)', Tri(A.A)) ≈ MyVector(b)' / Tri(A.A) + @test rdiv!(copy(B), Tri(A)) ≈ rdiv!(copy(B), Tri(A.A)) ≈ B / Tri(A.A) @test lmul!(Tri(A), copy(b)) ≈ lmul!(Tri(A.A), copy(b)) ≈ Tri(A.A) * MyVector(b) end @@ -154,29 +152,23 @@ BroadcastStyle(::Type{<:MyVector2{T}}) where {T} = MyBroadcastStyle{1}() @test cholesky(S, CRowMaximum()) \ b ≈ ldiv!(cholesky(Matrix(S), CRowMaximum()), copy(b)) @test cholesky(S) \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S)) \ b @test cholesky(S) \ b ≈ Symmetric(Matrix(S)) \ MyVector(b) - if VERSION >= v"1.9" - @test S \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S)) \ b - @test S \ b ≈ Symmetric(Matrix(S)) \ MyVector(b) - end + @test S \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S)) \ b + @test S \ b ≈ Symmetric(Matrix(S)) \ MyVector(b) S = Symmetric(MyMatrix(reshape(inv.(1:25),5,5) + 10I), :L) @test cholesky(S).U ≈ @inferred(cholesky!(deepcopy(S))).U @test cholesky(S,CRowMaximum()).U ≈ cholesky(Matrix(S),CRowMaximum()).U @test cholesky(S) \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S), :L) \ b @test cholesky(S) \ b ≈ Symmetric(Matrix(S), :L) \ MyVector(b) - if VERSION >= v"1.9" - @test S \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S), :L) \ b - @test S \ b ≈ Symmetric(Matrix(S), :L) \ MyVector(b) - end + @test S \ b ≈ Matrix(S) \ b ≈ Symmetric(Matrix(S), :L) \ b + @test S \ b ≈ Symmetric(Matrix(S), :L) \ MyVector(b) @testset "ldiv!" begin c = MyVector(randn(5)) @test ldiv!(lu(A), MyVector(copy(c))) ≈ A \ c @test_throws ErrorException ldiv!(eigen(randn(5,5)), c) @test ArrayLayouts.ldiv!(svd(A.A), Vector(c)) ≈ ArrayLayouts.ldiv!(similar(c), svd(A.A), c) ≈ A \ c - if VERSION ≥ v"1.8" - @test ArrayLayouts.ldiv!(similar(c), transpose(lu(A.A)), copy(c)) ≈ A'\c - end + @test ArrayLayouts.ldiv!(similar(c), transpose(lu(A.A)), copy(c)) ≈ A'\c B = Bidiagonal(randn(5), randn(4), :U) @test ldiv!(B, MyVector(copy(c))) ≈ B \ c @@ -254,24 +246,22 @@ BroadcastStyle(::Type{<:MyVector2{T}}) where {T} = MyBroadcastStyle{1}() @test_broken ldiv!(A, t) ≈ A\t @test ldiv!(A, copy(X)) ≈ A\X @test A\T ≈ A\T̃ - VERSION >= v"1.9" && @test A/T ≈ A/T̃ + @test 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 - VERSION >= v"1.9" && @test A/A ≈ I + @test A/A ≈ I @test A\MyVector(x) ≈ A\x @test A\MyMatrix(X) ≈ A\X - if VERSION >= v"1.9" - @test A/A ≈ A.A / A.A - @test x' / A ≈ x' / A.A - @test transpose(x) / A ≈ transpose(x) / A.A - @test transpose(x) / A isa Transpose - @test x' / A isa Adjoint - end + @test A/A ≈ A.A / A.A + @test x' / A ≈ x' / A.A + @test transpose(x) / A ≈ transpose(x) / A.A + @test transpose(x) / A isa Transpose + @test x' / A isa Adjoint @test D \ UpperTriangular(A) ≈ D \ UpperTriangular(A.A) @test UpperTriangular(A) \ D ≈ UpperTriangular(A.A) \ D @@ -349,9 +339,7 @@ BroadcastStyle(::Type{<:MyVector2{T}}) where {T} = MyBroadcastStyle{1}() B = randn(5,5) B̃ = MyMatrix(B) @test D*D ≈ Matrix(D)^2 - if VERSION ≥ v"1.7-" - @test D^2 ≈ D*D - end + @test D^2 ≈ D*D @test D*B ≈ Matrix(D)*B @test B*D ≈ B*Matrix(D) @test D*B̃ ≈ Matrix(D)*B̃ @@ -706,7 +694,7 @@ triangulardata(A::MyUpperTriangular) = triangulardata(A.A) @test_skip lmul!(U,view(copy(B),collect(1:5),1:5)) ≈ U * B @test MyMatrix(A) / U ≈ A / U - VERSION >= v"1.9" && @test U / MyMatrix(A) ≈ U / A + @test U / MyMatrix(A) ≈ U / A end @testset "* for infinite layouts" begin diff --git a/test/test_muladd.jl b/test/test_muladd.jl index e578468..dc6a239 100644 --- a/test/test_muladd.jl +++ b/test/test_muladd.jl @@ -104,9 +104,9 @@ Random.seed!(0) A, B = [1:4;], reshape([3:6;], 4, 1) D = Diagonal(Fill(3, 1)) M = MulAdd(2, A, D, 3, B) - @test copy(M) == (VERSION >= v"1.9" ? mul!(B, A, D, 2, 3) : 2 * A * D + 3 * B) + @test copy(M) == mul!(B, A, D, 2, 3) M = MulAdd(1, A, D, 0, B) - @test copy(M) == (VERSION >= v"1.9" ? mul!(B, A, D) : A * D) + @test copy(M) == mul!(B, A, D) end end @@ -674,7 +674,7 @@ Random.seed!(0) @test Q*Q' ≈ mul(Q,Q') @test Q'*Q' ≈ mul(Q',Q') @test Q'*Q ≈ mul(Q',Q) - VERSION >= v"1.8-" && @test Q*UpperTriangular(B) ≈ mul(Q, UpperTriangular(B)) + @test Q*UpperTriangular(B) ≈ mul(Q, UpperTriangular(B)) @test UpperTriangular(B)*Q' ≈ mul(UpperTriangular(B), Q') end @@ -838,19 +838,11 @@ Random.seed!(0) b = [quat(rand(4)...) for i in 1:4] c = [quat(rand(4)...) for i in 1:4, j in 1:1] M = MulAdd(α, b, D, β, c) - if VERSION >= v"1.9" - @test copy(M) ≈ mul!(copy(c), b, D, α, β) ≈ b * D * α + c * β - else - @test copy(M) ≈ b * D * α + c * β - end + @test copy(M) ≈ mul!(copy(c), b, D, α, β) ≈ b * D * α + c * β end @testset "Error paths" begin - if VERSION >= v"1.10.0" - Q = qr(rand(2,2), ColumnNorm()).Q - else - Q = qr(rand(2,2), Val(true)).Q - end + Q = qr(rand(2,2), ColumnNorm()).Q v = rand(Float32, 3) @test_throws DimensionMismatch ArrayLayouts.materialize!(ArrayLayouts.Rmul(v, Q)) @test_throws DimensionMismatch ArrayLayouts.materialize!(ArrayLayouts.Rmul(v, Q'))