first/isempty for RangeCumsum (#76)

* first/isempty for RangeCumsum * avoid ambiguities with Zeros multiplication * Revert "Revert "fill special cases"" This reverts commit 4482296. * increase cov * increase cov * tests pass * Update test_layoutarray.jl
JuliaLinearAlgebra · Sep 26, 2021 · 96bfa9d · 96bfa9d · dlfivefifty · Sep 26, 2021
1 parent 4482296
commit 96bfa9d
Show file tree

Hide file tree

Showing 7 changed files with 62 additions and 7 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,15 +1,15 @@
 name = "ArrayLayouts"
 uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
 authors = ["Sheehan Olver <solver@mac.com>"]
-version = "0.7.5"
+version = "0.7.6"
 
 [deps]
 FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 
 [compat]
-FillArrays = "0.11, 0.12"
+FillArrays = "0.12.6"
 julia = "1.5"
 
 [extras]

diff --git a/src/ArrayLayouts.jl b/src/ArrayLayouts.jl
@@ -162,6 +162,11 @@ getindex(A::LayoutVector, kr::Colon) = layout_getindex(A, kr)
 getindex(A::AdjOrTrans{<:Any,<:LayoutVector}, kr::Integer, jr::Colon) = layout_getindex(A, kr, jr)
 getindex(A::AdjOrTrans{<:Any,<:LayoutVector}, kr::Integer, jr::AbstractVector) = layout_getindex(A, kr, jr)
 
+*(a::Zeros{<:Any,2}, b::LayoutMatrix) = FillArrays.mult_zeros(a, b)
+*(a::LayoutMatrix, b::Zeros{<:Any,2}) = FillArrays.mult_zeros(a, b)
+*(a::Transpose{T, <:LayoutMatrix{T}} where T, b::Zeros{<:Any, 2}) = FillArrays.mult_zeros(a, b)
+*(a::Adjoint{T, <:LayoutMatrix{T}} where T, b::Zeros{<:Any, 2}) = FillArrays.mult_zeros(a, b)
+
 *(A::Diagonal{<:Any,<:LayoutVector}, B::Diagonal{<:Any,<:LayoutVector}) = mul(A, B)
 *(A::Diagonal{<:Any,<:LayoutVector}, B::AbstractMatrix) = mul(A, B)
 *(A::AbstractMatrix, B::Diagonal{<:Any,<:LayoutVector}) = mul(A, B)
@@ -306,6 +311,9 @@ end
 # printing
 ###
 
+const LayoutVecOrMat{T} = Union{LayoutVector{T},LayoutMatrix{T}}
+const LayoutVecOrMats{T} = Union{LayoutVecOrMat{T},SubArray{T,1,<:LayoutVecOrMat},SubArray{T,2,<:LayoutVecOrMat}}
+
 layout_replace_in_print_matrix(_, A, i, j, s) =
     i in colsupport(A,j) ? s : Base.replace_with_centered_mark(s)
 
@@ -315,7 +323,7 @@ Base.replace_in_print_matrix(A::Union{LayoutVector,
                                       UnitUpperTriangular{<:Any,<:LayoutMatrix},
                                       LowerTriangular{<:Any,<:LayoutMatrix},
                                       UnitLowerTriangular{<:Any,<:LayoutMatrix},
-                                      AdjOrTrans{<:Any,<:LayoutMatrix},
+                                      AdjOrTrans{<:Any,<:LayoutVecOrMat},
                                       HermOrSym{<:Any,<:LayoutMatrix},
                                       SubArray{<:Any,2,<:LayoutMatrix}}, i::Integer, j::Integer, s::AbstractString) =
     layout_replace_in_print_matrix(MemoryLayout(A), A, i, j, s)
@@ -339,8 +347,7 @@ include("cumsum.jl")
 # support overloading hcat/vcat for ∞-arrays
 ###
 
-const LayoutVecOrMat{T} = Union{LayoutVector{T},LayoutMatrix{T}}
-const LayoutVecOrMats{T} = Union{LayoutVecOrMat{T},SubArray{T,1,<:LayoutVecOrMat},SubArray{T,2,<:LayoutVecOrMat}}
+
 
 typed_hcat(::Type{T}, ::Tuple, A...) where T = Base._typed_hcat(T, A)
 typed_vcat(::Type{T}, ::Tuple, A...) where T = Base._typed_vcat(T, A)

diff --git a/src/cumsum.jl b/src/cumsum.jl
@@ -26,8 +26,10 @@ end
 
 Base.@propagate_inbounds getindex(c::RangeCumsum, kr::OneTo) = RangeCumsum(c.range[kr])
 
+first(r::RangeCumsum) = first(r.range)
 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))))

diff --git a/src/mul.jl b/src/mul.jl
@@ -182,6 +182,9 @@ macro layoutmul(Typ)
         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)
@@ -224,6 +227,7 @@ macro layoutmul(Typ)
             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)
 
             Base.:*(A::$Mod{<:Any,<:$Typ}, B::$Typ) = ArrayLayouts.mul(A,B)
             Base.:*(A::$Typ, B::$Mod{<:Any,<:$Typ}) = ArrayLayouts.mul(A,B)
@@ -267,6 +271,8 @@ dot(a, b) = materialize(Dot(a, b))
 @inline LinearAlgebra.dot(a::LayoutArray, b::LayoutArray) = dot(a,b)
 @inline LinearAlgebra.dot(a::LayoutArray, b::AbstractArray) = dot(a,b)
 @inline LinearAlgebra.dot(a::AbstractArray, b::LayoutArray) = dot(a,b)
+@inline LinearAlgebra.dot(a::LayoutVector, b::AbstractFill{<:Any,1}) = FillArrays._fill_dot(a,b)
+@inline LinearAlgebra.dot(a::AbstractFill{<:Any,1}, b::LayoutVector) = FillArrays._fill_dot(a,b)
 @inline LinearAlgebra.dot(a::LayoutArray{<:Number}, b::SparseArrays.SparseVectorUnion{<:Number}) = dot(a,b)
 @inline LinearAlgebra.dot(a::SparseArrays.SparseVectorUnion{<:Number}, b::LayoutArray{<:Number}) = dot(a,b)
 

diff --git a/src/muladd.jl b/src/muladd.jl
@@ -39,6 +39,8 @@ size(M::MulAdd) = map(length,axes(M))
 axes(M::MulAdd) = axes(M.C)
 
 similar(M::MulAdd, ::Type{T}, axes) where {T,N} = similar(Array{T}, axes)
+similar(M::MulAdd{<:DualLayout,<:Any,<:Any,<:Any,<:Adjoint}, ::Type{T}, axes) where {T,N} = similar(Array{T}, axes[2])'
+similar(M::MulAdd{<:DualLayout,<:Any,<:Any,<:Any,<:Transpose}, ::Type{T}, axes) where {T,N} = transpose(similar(Array{T}, axes[2]))
 similar(M::MulAdd, ::Type{T}) where T = similar(M, T, axes(M))
 similar(M::MulAdd) = similar(M, eltype(M))
 
@@ -371,6 +373,8 @@ scalarzero(::Type{<:AbstractArray{T}}) where T = scalarzero(T)
 
 fillzeros(::Type{T}, ax) where T<:Number = Zeros{T}(ax)
 mulzeros(::Type{T}, M) where T<:Number = fillzeros(T, axes(M))
+mulzeros(::Type{T}, M::Mul{<:DualLayout,<:Any,<:Adjoint}) where T<:Number = fillzeros(T, axes(M,2))'
+mulzeros(::Type{T}, M::Mul{<:DualLayout,<:Any,<:Transpose}) where T<:Number = transpose(fillzeros(T, axes(M,2)))
 
 # initiate array-valued MulAdd
 function _mulzeros!(dest::AbstractVector{T}, A, B) where T
@@ -412,4 +416,9 @@ function similar(M::MulAdd{<:Any,<:DualLayout,ZerosLayout}, ::Type{T}, (x,y)) wh
     @assert length(x) == 1
     trans = transtype(M.B)
     trans(similar(trans(M.B), T, y))
-end
+end
+
+const ZerosLayouts = Union{ZerosLayout,DualLayout{ZerosLayout}}
+copy(M::MulAdd{<:ZerosLayouts, <:ZerosLayouts, <:ZerosLayouts}) = M.C
+copy(M::MulAdd{<:ZerosLayouts, <:Any, <:ZerosLayouts}) = M.C
+copy(M::MulAdd{<:Any, <:ZerosLayouts, <:ZerosLayouts}) = M.C
diff --git a/test/test_cumsum.jl b/test/test_cumsum.jl
@@ -8,6 +8,7 @@ using ArrayLayouts, Test
         @test r[Base.OneTo(3)] == r[1:3]
         @test last(r) == r[end]
         @test diff(r) == diff(Vector(r))
+        @test first(r) == r[1]
     end
 
     a,b = RangeCumsum(Base.OneTo(5)), RangeCumsum(Base.OneTo(6))

diff --git a/test/test_layoutarray.jl b/test/test_layoutarray.jl
@@ -1,4 +1,4 @@
-using ArrayLayouts, LinearAlgebra, Test
+using ArrayLayouts, LinearAlgebra, FillArrays, Base64, Test
 import ArrayLayouts: sub_materialize
 
 if VERSION < v"1.7-"
@@ -173,6 +173,7 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor()
             A = MyMatrix(randn(5,5))
             b = randn(5)
             @test dot(b, A, b) ≈ b'*(A*b) ≈ b'A*b
+            @test dot(b, A, b) ≈ transpose(b)*(A*b) ≈ transpose(b)A*b
         end
 
         @testset "dual vector * symmetric (#40)" begin
@@ -265,8 +266,32 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor()
         @test Transpose(T) * A ≈ Transpose(T) * A.A
         @test Transpose(T)A' ≈ Adjoint(T)A' ≈ Adjoint(T)Transpose(A) ≈ Transpose(T)Transpose(A)
         @test Transpose(A)Adjoint(T) ≈ A'Adjoint(T) ≈ A'Transpose(T) ≈ Transpose(A)Transpose(T)
+
+        @test Zeros(5)' * A ≡ Zeros(5)'
+        @test transpose(Zeros(5)) * A ≡ transpose(Zeros(5))
+
+        @test A' * Zeros(5) ≡ Zeros(5)
+        @test Zeros(3,5) * A ≡ Zeros(3,5)
+        @test A * Zeros(5,3) ≡ Zeros(5,3)
+        @test A' * Zeros(5,3) ≡ Zeros(5,3)
+        @test transpose(A) * Zeros(5,3) ≡ Zeros(5,3)
+        @test A' * Zeros(5) ≡ Zeros(5)
+        @test transpose(A) * Zeros(5) ≡ Zeros(5)
+
+        b = MyVector(randn(5))
+        @test A' * b ≈ A' * b.A
     end
 
+    @testset "AbstractQ" begin
+        A = MyMatrix(randn(5,5))
+        Q = qr(randn(5,5)).Q
+        @test Q'*A ≈ Q'*A.A
+        @test Q*A ≈ Q*A.A
+        @test A*Q ≈ A.A*Q
+        @test A*Q' ≈ A.A*Q'
+    end
+
+
     @testset "concat" begin
         a = MyVector(randn(5))
         A = MyMatrix(randn(5,6))
@@ -278,6 +303,11 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor()
         @test [a A A] == [Array(a) Array(A) Array(A)]
         @test [a Array(A) A] == [Array(a) Array(A) Array(A)]
     end
+
+    @testset "dot" begin
+        a = MyVector(randn(5))
+        @test dot(a, Zeros(5)) ≡ dot(Zeros(5), a) ≡ 0.0
+    end
 end
 
 struct MyUpperTriangular{T} <: AbstractMatrix{T}