Fix Vector'Diagonal (and Transpose as well) to avoid infinite recursion. #26924
Conversation
andreasnoack
changed the title
stdlib/LinearAlgebra/src/diagonal.jl
Outdated
| @@ -481,6 +481,10 @@ end | |||
| # dismabiguation methods: * of Diagonal and Adj/Trans AbsVec | |||
| *(A::Diagonal, B::Adjoint{<:Any,<:AbstractVector}) = A * copy(B) | |||
| *(A::Diagonal, B::Transpose{<:Any,<:AbstractVector}) = A * copy(B) | |||
| *(A::Adjoint{<:Any,<:AbstractVector}, B::Diagonal) = copy(A) * B | |||
| *(A::Transpose{<:Any,<:AbstractVector}, B::Diagonal) = copy(A) * B | |||
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal) = Adjoint(map((t,s) -> t'*s, D.diag, parent(x))) | |||
There was a problem hiding this comment.
Perhaps (D' * x')' would be simpler / clearer, if methods appropriate for D' * x' exist? (Edit: I.e. methods for *(Adjoint{<:Any,<:Diagonal}, AbstractVector).)
stdlib/LinearAlgebra/src/diagonal.jl
Outdated
| *(A::Transpose{<:Any,<:AbstractVector}, B::Diagonal) = copy(A) * B | ||
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal) = Adjoint(map((t,s) -> t'*s, D.diag, parent(x))) | ||
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal, y::AbstractVector) = | ||
| mapreduce(t -> t[1]'*t[2]*t[3], +, zip(parent(x), D.diag, y)) |
There was a problem hiding this comment.
Alternatively mapreduce(t -> t[1]*t[2]*t[3], +, zip(x, D.diag, y))?
stdlib/LinearAlgebra/src/diagonal.jl
Outdated
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal) = Adjoint(map((t,s) -> t'*s, D.diag, parent(x))) | ||
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal, y::AbstractVector) = | ||
| mapreduce(t -> t[1]'*t[2]*t[3], +, zip(parent(x), D.diag, y)) | ||
| *(A::Transpose{<:Any,<:AbstractVector}, B::Diagonal) = Transpose(map(*, D.diag, parent(x))) |
There was a problem hiding this comment.
B in this method signature should be D? (Perhaps test?)
Similarly to the above, perhaps Transpose(transpose(D) * transpose(x)) would be simpler / clearer, if methods appropriate for transpose(D) * transpose(x) exist? (Edit: I.e. methods for *(Transpose{<:Any,<:Diagonal}, AbstractVector).)
stdlib/LinearAlgebra/src/diagonal.jl
Outdated
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal, y::AbstractVector) = | ||
| mapreduce(t -> t[1]'*t[2]*t[3], +, zip(parent(x), D.diag, y)) | ||
| *(A::Transpose{<:Any,<:AbstractVector}, B::Diagonal) = Transpose(map(*, D.diag, parent(x))) | ||
| *(x::Adjoint{<:Any,<:AbstractVector}, D::Diagonal, y::AbstractVector) = |
There was a problem hiding this comment.
This method signature is identical to the signature above? Perhaps you meant *(x::Transpose{<:Any,<:AbstractVector}, D::Diagonal, y::AbstractVector) instead? In such case the method body seems slightly amiss as well, particularly for nested objects? Perhaps mapreduce(t -> t[1]*t[2]*t[3], +, zip(x, D.diag, y)) instead?
ararslan
changed the title
Also add optimized methods for x'D*y to avoid allocating temporary vector
andreasnoack
changed the title
|
Looks like this is working now. |
Also add optimized methods for x'D*y to avoid allocating temporary vector
Fixes JuliaLang/LinearAlgebra.jl#517