Skip to content

generic_matvecmul fails for some non-scalar element type and adjoint combinations #45029

Open

Description

Here is a super stripped down MWE.

In the example below, I would expect M'*x to return a Vector{Float64} equivalent to [2.0]. But it looks like the fact that the matrix elements themselves should be adjoint is lost somewhere along the way to generic_matvecmul, which attempts straight multiplication of two vectors, rather than taking their inner product.

julia> M = [[[1.0]];;]
1×1 Matrix{Vector{Float64}}:
 [1.0]

julia> x = [[[2.0]];]
1-element Vector{Vector{Float64}}:
 [2.0]

julia> M'*x
ERROR: MethodError: no method matching *(::Vector{Float64}, ::Vector{Float64})
Closest candidates are:
  *(::Any, ::Any, ::Any, ::Any...) at /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/base/operators.jl:655
  *(::StridedMatrix{T}, ::StridedVector{S}) where {T<:Union{Float32, Float64, ComplexF32, ComplexF64}, S<:Real} at /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:44
  *(::StridedVecOrMat, ::Adjoint{<:Any, <:LinearAlgebra.LQPackedQ}) at /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/lq.jl:266
  ...
Stacktrace:
 [1] generic_matvecmul!(C::Vector{Union{}}, tA::Char, A::Matrix{Vector{Float64}}, B::Vector{Vector{Float64}}, _add::LinearAlgebra.MulAddMul{true, true, Bool, Bool})
   @ LinearAlgebra /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:745
 [2] mul!
   @ /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:129 [inlined]
 [3] mul!
   @ /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:275 [inlined]
 [4] *(adjA::Adjoint{Adjoint{Float64, Vector{Float64}}, Matrix{Vector{Float64}}}, x::Vector{Vector{Float64}})
   @ LinearAlgebra /Applications/Julia-1.7.app/Contents/Resources/julia/share/julia/stdlib/v1.7/LinearAlgebra/src/matmul.jl:113
 [5] top-level scope
   @ REPL[152]:1

In contrast, these three other similar operations, seem to work correctly:

julia> M'*M
1×1 Matrix{Float64}:
 1.0

julia> x'*x
4.0

julia> N = [[[1.0]'];;]
1×1 Matrix{Adjoint{Float64, Vector{Float64}}}:
 [1.0;;]

julia> N*x
1-element Vector{Float64}:
 2.0
julia> versioninfo()
Julia Version 1.7.2
Commit bf53498635 (2022-02-06 15:21 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin19.5.0)
  CPU: Intel(R) Core(TM) i7-9750H CPU @ 2.60GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-12.0.1 (ORCJIT, skylake)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions