fix matrix multiplication interaction with HermOrSym and Diagonal (#2…

…2474)

(cherry picked from commit 0990d9b)

base/linalg/diagonal.jl

@@ -150,9 +150,9 @@ end
 (*)(D::Diagonal, B::AbstractTriangular) = A_mul_B!(D, copy(B))
 
 (*)(A::AbstractMatrix, D::Diagonal) =
-    scale!(similar(A, promote_op(*, eltype(A), eltype(D.diag))), A, D.diag)
+    scale!(similar(A, promote_op(*, eltype(A), eltype(D.diag)), size(A)), A, D.diag)
 (*)(D::Diagonal, A::AbstractMatrix) =
-    scale!(similar(A, promote_op(*, eltype(A), eltype(D.diag))), D.diag, A)
+    scale!(similar(A, promote_op(*, eltype(A), eltype(D.diag)), size(A)), D.diag, A)
 
 A_mul_B!(A::Union{LowerTriangular,UpperTriangular}, D::Diagonal) =
     typeof(A)(A_mul_B!(A.data, D))

test/linalg/diagonal.jl

@@ -346,3 +346,18 @@ end
     @test logm(D) == Diagonal([logm([1 2; 3 4]), logm([1 2; 3 4])])
     @test sqrtm(D) == Diagonal([sqrtm([1 2; 3 4]), sqrtm([1 2; 3 4])])
 end
+
+@testset "multiplication with Symmetric/Hermitian" begin
+    for T in (Float64, Complex128)
+        D = Diagonal(randn(T, n))
+        A = randn(T, n, n); A = A'A
+        S = Symmetric(A)
+        H = Hermitian(A)
+        for f in (*, Ac_mul_B, A_mul_Bc, Ac_mul_Bc, At_mul_B, A_mul_Bt, At_mul_Bt)
+            @test f(D, S) ≈ f(Matrix(D), Matrix(S))
+            @test f(D, H) ≈ f(Matrix(D), Matrix(H))
+            @test f(S, D) ≈ f(Matrix(S), Matrix(D))
+            @test f(S, H) ≈ f(Matrix(S), Matrix(H))
+        end
+    end
+end