More mul! overloads needed by BandedBlockBandedMatrix (#165)

* Update mul.jl

* systematic overloading

* Update test_layoutarray.jl

* ambiguities with Fill multiplication

* Update test_layoutarray.jl

* v0.2.1

* Update mul.jl

* Update test_layoutarray.jl

Author: dlfivefifty

Project.toml

@@ -1,7 +1,7 @@
 name = "ArrayLayouts"
 uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
 authors = ["Sheehan Olver <[email protected]>"]
-version = "1.2.0"
+version = "1.2.1"
 
 [deps]
 FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"

src/ArrayLayouts.jl

@@ -203,6 +203,14 @@ getindex(A::AdjOrTrans{<:Any,<:LayoutVector}, kr::Integer, jr::AbstractVector) =
 *(a::Adjoint{T, <:LayoutMatrix{T}} where T, b::Zeros{<:Any, 2}) = FillArrays.mult_zeros(a, b)
 *(A::Adjoint{<:Any, <:Zeros{<:Any,1}}, B::Diagonal{<:Any,<:LayoutVector}) = (B' * A')'
 *(A::Transpose{<:Any, <:Zeros{<:Any,1}}, B::Diagonal{<:Any,<:LayoutVector}) = transpose(transpose(B) * transpose(A))
+*(a::Adjoint{<:Number,<:LayoutVector}, b::Zeros{<:Number,1})= FillArrays._adjvec_mul_zeros(a, b)
+function *(a::Transpose{T, <:LayoutVector{T}}, b::Zeros{T, 1}) where T<:Real
+    la, lb = length(a), length(b)
+    if la ≠ lb
+        throw(DimensionMismatch("dot product arguments have lengths $la and $lb"))
+    end
+    return zero(T)
+end
 
 *(A::Diagonal{<:Any,<:LayoutVector}, B::Diagonal{<:Any,<:LayoutVector}) = mul(A, B)
 *(A::Diagonal{<:Any,<:LayoutVector}, B::AbstractMatrix) = mul(A, B)

src/mul.jl

@@ -317,10 +317,20 @@ end
 
 
 # mul! for subarray of layout matrix
-LinearAlgebra.mul!(C::SubArray{<:Any,2,<:LayoutMatrix}, A::SubArray{<:Any,2,<:LayoutMatrix}, B::SubArray{<:Any,2,<:LayoutMatrix}, α::Number, β::Number) =
+const SubLayoutMatrix = Union{SubArray{<:Any,2,<:LayoutMatrix}, SubArray{<:Any,2,<:AdjOrTrans{<:Any,<:LayoutMatrix}}}
+
+LinearAlgebra.mul!(C::AbstractMatrix, A::SubLayoutMatrix, B::AbstractMatrix, α::Number, β::Number) =
+    ArrayLayouts.mul!(C, A, B, α, β)
+LinearAlgebra.mul!(C::AbstractMatrix, A::AbstractMatrix, B::SubLayoutMatrix, α::Number, β::Number) =
+    ArrayLayouts.mul!(C, A, B, α, β)
+LinearAlgebra.mul!(C::AbstractMatrix, A::SubLayoutMatrix, B::LayoutMatrix, α::Number, β::Number) =
     ArrayLayouts.mul!(C, A, B, α, β)
-LinearAlgebra.mul!(C::AbstractMatrix, A::SubArray{<:Any,2,<:LayoutMatrix}, B::SubArray{<:Any,2,<:LayoutMatrix}, α::Number, β::Number) =
+LinearAlgebra.mul!(C::AbstractMatrix, A::LayoutMatrix, B::SubLayoutMatrix, α::Number, β::Number) =
     ArrayLayouts.mul!(C, A, B, α, β)    
+LinearAlgebra.mul!(C::AbstractMatrix, A::SubLayoutMatrix, B::SubLayoutMatrix, α::Number, β::Number) =
+    ArrayLayouts.mul!(C, A, B, α, β)
+LinearAlgebra.mul!(C::AbstractVector, A::SubLayoutMatrix, B::AbstractVector, α::Number, β::Number) =
+    ArrayLayouts.mul!(C, A, B, α, β)
 
 
 ###

src/muladd.jl

@@ -224,9 +224,9 @@ function _default_blasmul!(::IndexCartesian, α, A::AbstractMatrix, B::AbstractV
     z = zero(A[1,1]*B[1] + A[1,1]*B[1])
 
     @inbounds for k in colsupport(B,1)
-        b = B[k]
+        b = α * B[k]
         for i = colsupport(A,k)
-            C[i] += α * A[i,k] * b
+            C[i] += A[i,k] * b
         end
     end
 

test/test_layoutarray.jl

@@ -361,6 +361,11 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor()
 
         b = MyVector(randn(5))
         @test A' * b ≈ A' * b.A
+
+        @test b'*Zeros(5) == 0
+        @test transpose(b)*Zeros(5) == 0
+        @test_throws DimensionMismatch b'*Zeros(6)
+        @test_throws DimensionMismatch transpose(b)*Zeros(6)
     end
 
     @testset "AbstractQ" begin
@@ -474,6 +479,28 @@ MemoryLayout(::Type{MyVector}) = DenseColumnMajor()
             @test LinearAlgebra.copymutable_oftype(A, BigFloat) == A
         end
     end
+
+    @testset "mul! with subarrays" begin
+            A = MyMatrix(randn(3,3))
+            V = view(A, 1:3, 1:3)
+            B = randn(3,3)
+            x = randn(3)
+            @test mul!(similar(B), V, B) ≈ A * B
+            @test mul!(similar(B), B, V) ≈ B * A
+            @test mul!(similar(B), V, V) ≈ A^2
+            @test mul!(similar(B), V, A) ≈ A * A
+            @test mul!(similar(B), A, V) ≈ A * A
+            @test mul!(MyMatrix(randn(3,3)), A, V) ≈ A * A
+            @test mul!(similar(x), V, x) ≈ V * x
+
+            @test mul!(copy(B), V, B, 2.0, 3.0) ≈ 2A * B + 3B
+            @test mul!(copy(B), B, V, 2.0, 3.0) ≈ 2B * A + 3B
+            @test mul!(copy(B), V, V, 2.0, 3.0) ≈ 2A^2 + 3B
+            @test mul!(copy(B), V, A, 2.0, 3.0) ≈ 2A * A + 3B
+            @test mul!(copy(B), A, V, 2.0, 3.0) ≈ 2A * A + 3B
+            @test mul!(MyMatrix(copy(B)), A, V, 2.0, 3.0) ≈ 2A * A + 3B
+            @test mul!(copy(x), V, x, 2.0, 3.0) ≈ 2A * x + 3x
+    end
 end
 
 struct MyUpperTriangular{T} <: AbstractMatrix{T}