Request for `inv` and `\` for `Diagonal`

[jeksterslab]

Requesting for symbolic inv and \ to be available for Diagonal matrix. The example below returns errors for inv and \ for operations involving a Diagonal matrix D.

julia> using LinearAlgebra

julia> using Symbolics

julia> @variables a b c
(a, b, c)

julia> A = [
         0 a 1;
         0 0 0;
         0 0 0;
       ]
3×3 Array{Num,2}:
 0  a  1
 0  0  0
 0  0  0

julia> B = [
         0 0 0;
         0 b 0;
         0 0 c;
       ]
3×3 Array{Num,2}:
 0  0  0
 0  b  0
 0  0  c

julia> C = (I-A)\B/(I-A)'
3×3 Array{Num,2}:
 c + b*(a^2)  a*b    c
         a*b    b    0.0
           c    0.0  c

julia> D = Diagonal(sqrt.(diag(C)))
3×3 Diagonal{Num,Array{Num,1}}:
 sqrt(c + b*(a^2))        ⋅        ⋅
                 ⋅  sqrt(b)        ⋅
                 ⋅        ⋅  sqrt(c)

julia> inv(D)
ERROR: TypeError: non-boolean (Num) used in boolean context
Stacktrace:
 [1] inv(::Diagonal{Num,Array{Num,1}}) at /build/julia/src/julia-1.5.3/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/diagonal.jl:606
 [2] top-level scope at REPL[106]:1
 [3] run_repl(::REPL.AbstractREPL, ::Any) at /build/julia/src/julia-1.5.3/usr/share/julia/stdlib/v1.5/REPL/src/REPL.jl:288

julia> D\C/D
ERROR: TypeError: non-boolean (Num) used in boolean context
Stacktrace:
 [1] ldiv!(::Diagonal{Num,Array{Num,1}}, ::Array{Num,2}) at /build/julia/src/julia-1.5.3/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/diagonal.jl:589
 [2] \(::Diagonal{Num,Array{Num,1}}, ::Array{Num,2}) at /build/julia/src/julia-1.5.3/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/diagonal.jl:597
 [3] top-level scope at REPL[107]:1
 [4] run_repl(::REPL.AbstractREPL, ::Any) at /build/julia/src/julia-1.5.3/usr/share/julia/stdlib/v1.5/REPL/src/REPL.jl:288

┆Issue is synchronized with this Trello card by Unito

[ChrisRackauckas]

I think it just needs an overload, and the overload there is pretty trivial since it's just elementwise. inv(D::Diagonal{Num}) = inv.(D), and similarly for the ldiv! overload.

[jeksterslab]

I am not sure how to overload a function properly but I tried the following for inv. I am still figuring out something for \. I think some simplification rules can be added to clean-up the final result.

julia> function inv(D::Diagonal{Num})
       Diagonal(Base.inv.(diag(D)))
       end
inv (generic function with 1 method)

julia> inv(D)
3×3 Diagonal{Num,Array{Num,1}}:
 (sqrt(c + b*(a^2)))^-1             ⋅             ⋅
                      ⋅  (sqrt(b))^-1             ⋅
                      ⋅             ⋅  (sqrt(c))^-1

julia> inv(D) * C * inv(D)'
3×3 Array{Num,2}:
    (c + (b*(a^2)))*((sqrt(c + b*(a^2)))^-2)  a*b*((sqrt(b))^-1)*((sqrt(c + b*(a^2)))^-1)    c*((sqrt(c))^-1)*((sqrt(c + b*(a^2)))^-1)
 a*b*((sqrt(b))^-1)*((sqrt(c + b*(a^2)))^-1)                             b*((sqrt(b))^-2)                                            0.0
   c*((sqrt(c))^-1)*((sqrt(c + b*(a^2)))^-1)                                            0.0                           c*((sqrt(c))^-2)

[jeksterslab]

I tried the following for ldiv! but I am not sure if the function definition is appropriate.

julia> function ldiv!(D::Diagonal{Num}, B)
         inv(D) * B
       end
ldiv! (generic function with 1 method)

julia> ldiv!(D, C) / D
3×3 Array{Num,2}:
    (c + (b*(a^2)))*((sqrt(c + b*(a^2)))^-2)  a*b*((sqrt(b))^-1)*((sqrt(c + b*(a^2)))^-1)    c*((sqrt(c))^-1)*((sqrt(c + b*(a^2)))^-1)
 a*b*((sqrt(b))^-1)*((sqrt(c + b*(a^2)))^-1)                             b*((sqrt(b))^-2)                                            0.0
   c*((sqrt(c))^-1)*((sqrt(c + b*(a^2)))^-1)                                            0.0                           c*((sqrt(c))^-2)

[ChrisRackauckas]

You might get a simpler representation by directly doing lu(D) \ B.