Symbolics.jl, few reflections of beta user

Part of the March I spend learning basics of Symbolics.jl package and throwing questions and commits with one line changed (I now know, you shouldn’t do that) on JuliaSymbolics team. I want to share few reflections about it, hoping that they would be starting point for discussions.

Next two paragraph will be quit long and probably boring, yet without knowing why I want to achieve, my point of view will be hard to understand. I simplified some technical stuff and omitted few things, still they are quit long.

One of my pet projects is to assemble framework for solving numerically equation of motion of point particle on given space-time in the settings of general relativity theory. DifferentialEquations.jl give me already all I want at this moment and much more. But as anyone who study this subject know, even if we pull out from our magical hat reasonable form of metric tensor, which describes effects of gravitation on space-time, writing down equations of motion for such can be long and error prone procedure. Since changing coordinate system is such context is normal things to do, for many reasons, symbolic computations is what you often want to have at hand. Heaving that as motivation I checked Symbolics.jl and SymbolicTensors.jl. The latter was written with general relativity in mind, so this is just what I want. Cheers for Robbie Rosati and @miguelraz. In this moment the best thing I can do is to try learn something SymbolicUtils.jl and try to switch SymbolicTensors.jl backend from SymPy to SymbolicUtils.jl and benchmark performance. After that, we will see what to do.

In result I want to write few reflections on Symbolics.jl from the point of view of user that go through tutorial. We have here true Computer Algebra Systems (CAS) expersts like @rfateman, @carette and @dimpase, so I don’t even dear to pretend that I know about this subject more than I learn reading recent discussions on Julia Discourses and Julia. I just want to write down what I found as the user and few reflections that come to my mind along the way. They’re probably quit shallow, but I hope someone much better will correct them and we can learn from his/her knowledge.

  1. Symbolics.jl is very beta. Many things works in odd way or don’t work at all. But, JuliaSymbolics team is aware of this and works hard to make it better.

  2. Central concept of Symbolics.jl seems to be symbolic variable: @variables x. This concept looks similar to dependent and independent variables (inconvenient names in my opinion) and some of this properties supports that claim.

@variables x y

z = x^2 + y
typeof(z)

result in Num. According to tutorial "@variables returns Sym or istree objects wrapped in Num in order to make them behave like subtypes of Real" and documentation states that “Sym, Term, and FnType are from SymbolicUtils.jl. Note that in Symbolics, we always use Sym{Real}, Term{Real}, and FnType{Tuple{Any}, Real}.” I’m already a little bit lost.

I don’t have good understanding notion of symbols in Julia, but I have some materials to read about this topic. Second, Julia have Lisp-like symbols :x and typeof(:x) returns Symbol. On the other hand SymbolicUtils.jl it have it own notion of symbols defined using macro @syms. This make me confused about precise meaning of terms “symbol” and “symbolic variable” in this context.

  1. Consider that we now add code
A = [x   0 y
     y   x 0
     1.0 x y^2]

A have now type Matrix{Num} (alias for Array{Num, 2}). This results are make me feel uneasy.

From what I find Num(val) is a part of Symbolics.jl, its description is Wrap anything in a type that is a subtype of Real. According to documentation of Julia we have in core language Real <: Number, real numbers are subtype of numbers, I have no problem with that. But this suggest to me that numbers ARE subtypes or real numbers, with which I disagree. My opinion is that p-adic numbers are as good numbers as real numbers, since they are both completions of the rational numbers, but p-adic numbers are NOT a subtype of real numbers. This probably isn’t important for general relativity, but just raises objections.

Second thing that bother me is that this choice force symbolic variable to represent real number, which can be to limiting. If I understand correctly what @dimpase said previously, many practical calculations involves p-adic numbers in the intermediate steps. And thought that in the future need of such algorithms can lead to implementing p-adic numbers as subtypes of Real make me cringe.

Also, if symbolic variable is Symbolics.jl notion of symbol in algebraic system, why there is a division between symbolic variable and matrix of symbolic variable? My philosophical opinion, so it is very questionable thing, is that symbols should denote whole range of mathematical objects (to consider only such class of things) like numbers, matrices, functions, sets, and so on. Also, in physics manipulation of such “composed objects” as simple symbols is omnipresent, so this make me worried is this a proper way to go. Maybe there are important theoretical and technical reasons to do this in such, I’m very ignorant about CAS, at this moment I feel big uneasy about this.

  1. As @carette says on Zulip “symbols are not imperative variable”. Consider this code.
@variables x y

z = x^2 + y
# `x`, `y` and `z` have type `Num`

x = 1 # `x` have type `Int64`
z     # have value (?) `x^2 - y` and type `Num`
z - x # have type `Num` and value `y + x^2 - 1`

This fortunately don’t suffer from problem that @carette mention, which is that in Maple after putting x = 1 value of z is changed to y + 1, but this look to me as loophole for bugs. Again, I maybe just to ignorant to understand what this works this way.

You can make desirable substitution by

substitute(x^3 + y^2 + z/2, Dict(x => 3, y => 2, z => 1))

Which works very fine for me.

  1. Interplay between symbolic computation and standard code can works very smooth, but current beta version have many rough corners here. But ability to define function f(x) in Julia and arrive at it symbolic forms by just calling it with symbolic variable is a great fun.

  2. Turning symbolic code after some symbolic reorganization into compiled numerical code what I dream of for point particle on given space-time. For it purpose symbols variables being subtypes of reals are probably enough (but who can be sure, physics is strange creature), so if we can improve SymbolicsTensors.jl, my aim can be not so far, far away. Again some rough edges are present

I probably should write more my reflections, but this is already quite long post and on Zulip there are new, very important threads about CAS, which probably make me rething few things, so I post it in this form. I hope that someone learn something from my experience. If I find good reasons I will write more on this topic.

3 Likes

Thanks! Most of these questions are actually related to the WIP PRs.

Num is just a hack for a compiler hook and not intended to be a lasting choice. The real solution is to get Sym{T} <: T or some symbolic tracer for such dispatch control:

What’s going on is that a metadata system was just implemented, so Sym does have a type, Sym{T}, with other information. But Sym{T} then broadcasts like Sym <: Any. So if you call a generic Julia function that wants a number, then f(x) fails. Num <: Real and so Num(x::Sym) <: Real and thus f(x) works by default, and it traces through a bunch of code. That gets us going to day, but what we really need is for f(x::Sym{T}) to note T <: Real and take the right dispatch. The proposed syntax is then that:

@variables x::Complex y z[1:2,1:2]::Quaternion

generates the expected stuff, with the default just being that it’s a symbolic real (with the behavior currently matching Num. The difficulty is that we need some compiler hooks that don’t exist, so Num is a hack we probably have to live with until Julia v1.7. However, I know Keno is running into this same issue with Diffractor.jl, so we’re brainstorming some real solutions.

So yes, your worries about things being real is simply because you haven’t seen the fact that you can change the type (but it just doesn’t work right now :wink:).

That’s made to solve this issue. Making array of structs was just a short term solution to get the downstream tests going.

Symbols and imperative variables purposefully are different, as otherwise metaprogramming would be impossible without eval. The tl;dr of @variables x is that it creates a Julia variable x with a Num(Sym(:x)). User programs use macros so the two are always named the same, but there’s purposefully a separation because of programs that act on symbolic variables. For example, you want to say x = variable_from_user, and still have x refer to :z if what the user passed into the function was defined by @variables z. Otherwise you’d always have to match global names. And that makes it clear why it solves the problem of Maple (and why writing Maple would be ridiculously hard because of violating this!).

8 Likes

Thank you for answering my worries. This post would be better if I dig more into Symbolics.jl, but as I explained what my pet project needs now is SymbolicUtils.jl, so I will be messing around (hopeful not literally :wink:) in that part of ecosystem.

Also @carette is posting now his reflections on Zulip, which are streatching to maximum my \lambda calculus and algebra knowledge, so you have true CAS expertise to evaluate your project choices.

Sound very exiting. I wish you luck.

3 Likes