Cannot substitute complex value into complex symbolic expression

[akriegman]

Minimal working example:

julia> using Symbolics

julia> @variables z
1-element Vector{Num}:
 z

julia> substitute(im*z, z => im)
ERROR: InexactError: Real(im)

When I would have expected the result to be -1. I took a dive into the code and figured out what's happening. im*z is represented as Complex{Num}(0, z). When we try to substitute im for z we make a call to Complex{Real}(re::Bool, im::Complex{Bool}), since that's the most literal way to replace z with im. This throws an error since complex numbers in Julia have to have real coefficients.

One potential fix could be to represent im*z as Mul(im, z) instead. Then when we substitute a complex value for z we would multiply that value instead of trying to make it the coefficient of a complex number. I assume this isn't the easiest solution though.

Maybe we could have symbolic variables that are considered to be complex. Then im*x where x is a normal Sym would be Complex(0, x) like how it works now, while im*z where z is <: Complex would be Mul(im, z), since Complex(0, z) would be an error, as z isn't <: Real in the first place. Or maybe Symbolics already has a better way of dealing with complex numbers that I don't know about.

┆Issue is synchronized with this Trello card by Unito

[akriegman]

Update: I was able to get around this by defining the following overloads:

Base.Complex{T}(r::Number, i::Complex) where {T<:Real} = r + im * i
Base.Complex{T}(r::Complex, i::Number) where {T<:Real} = r + im * i
Base.Complex{T}(r::Complex, i::Complex) where {T<:Real} = r + im * i

Now everything works as expected. This is a little bit hacky, but is also arguably what Complex should do when given a complex argument anyways.

[michikawa07]

@akriegman Do you have a pull request? Or do you plan to do so?

[adomasbaliuka]

This might be related to an issue where substituting variables inside a complex matrix does not work at all (has no effect and leaves the variables as they are).

julia> @variables a
1-element Vector{Num}:
 a

julia> M = [im*a 0
       a 0]
2×2 Matrix{Complex{Num}}:
 a*im  0
    a  0

julia> substitute(M, Dict([a=>1.]))
2×2 Matrix{Complex{Num}}:
 a*im  0
    a  0