search

JuliaSymbolics / Symbolics.jl

Symbolic programming for the next generation of numerical software
https://docs.sciml.ai/Symbolics/stable/
Other
1.37k stars 153 forks source link

Eigen and eigvals #259

Open baggepinnen opened 3 years ago

baggepinnen commented 3 years ago

SymPy handles eigenvalues of symbolic matrices, provided that they are not too large. Loading GenericLinearAlgebra and trying the same with Symbolics results in 

ERROR: TypeError: non-boolean (Num) used in boolean context
julia> using Symbolics, GenericLinearAlgebra

julia> @variables Jm J0

2-element Vector{Num}:
 Jm
 J0

julia> eigvals([Jm J0; 0 1])
ERROR: TypeError: non-boolean (Num) used in boolean context
Stacktrace:
 [1] ishermitian(A::Matrix{Num})
   @ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1248
 [2] eigvals!(A::Matrix{Num}; sortby::typeof(LinearAlgebra.eigsortby), kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ GenericLinearAlgebra ~/.julia/packages/GenericLinearAlgebra/z4KPO/src/eigenGeneral.jl:257
 [3] eigvals!
   @ ~/.julia/packages/GenericLinearAlgebra/z4KPO/src/eigenGeneral.jl:257 [inlined]
 [4] eigvals(A::Matrix{Num}; kws::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/eigen.jl:326
 [5] eigvals(A::Matrix{Num})
   @ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/eigen.jl:326

julia> using SymPy

julia> @vars Jm J0
(Jm, J0)

julia> eigvals([Jm J0; 0 1])
2-element Vector{Sym}:
 Jm
  1

┆Issue is synchronized with this Trello card by Unito

YingboMa commented 3 years ago

Basically, dup of #220

YingboMa commented 3 years ago

That's not how GenericLinearAlgebra works BTW...

baggepinnen commented 3 years ago

That's not how GenericLinearAlgebra works BTW...

Not sure what you mean by that, GenericLinearAlgebra does a have a method of eigvals which is generic wrt the matrix element type, so it felt natural to try that, granted that it's most likely targeted at numerical computing only.

ChrisRackauckas commented 3 years ago

The generic algorithm isn't going to work. What makes this special is that eigenvalues are solutions to the characteristic polynomials, and so you know by Abel's theorem that you definitely cannot do this for NxN with N>=5 as if you could, then you'd have an algorithm for solving higher order polynomials effectively. Something could be specialized to small matrices, but you'd probably what to do it via the polynomial form.

https://en.wikipedia.org/wiki/Eigenvalue_algorithm#Algorithms

Numerical eigenvalue solvers use an algorithm that iterates to tolerance, but you just can't do that if you have symbolic variables because symbolic control flow cannot be evaluated without values. The error:

ERROR: TypeError: non-boolean (Num) used in boolean context

shows up when you're doing something like if x < tol. If x is a number, then this is true or false. If x is symbol, then it's x < tol, so Julia just cannot know how many iterations to do and throws an error.

Hopefully that makes the error, the reason for the error, the reason why it has to be there in general, and exactly what the workaround would look like (i.e. hardcoding polynomial solves for up to the quartic and using that). But the moral of the story is, if you see that error, you're likely doing a non-quasi-static algorithm, i.e. some algorithm which requires knowing the numerical values to know the number of steps, and those algorithms are by their very nature not compatible with symbolic computing.