# CAS Best Practices

**URL:** <https://discourse.julialang.org/t/cas-best-practices/58092>\
**Category:** Machine Learning\
**Tags:** symbolic\
**Created:** [March 27, 2021, 8:51pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092 "2021-03-27T20:51:11Z")\
**Posts on this page:** 20\
**Page:** 4

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [April 4, 2021, 3:08pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/64 "2021-04-04T15:08:28Z")

</div>

> [@ChrisRackauckas](#):
>
> It’s zero(x) or zero(T) in Julia

This doesn’t work if `x` is an array — `zero(x) * x` will throw if `x` is not a square matrix. In the context of rings it is fine, of course, though it still won’t be equivalent to `x-x` for square matrices that may contain non-finite values.

I guess you could do something like `zeroscalar(x) * x` where

```julia
zeroscalar(x) = zero(x)
zeroscalar(::Type{<:AbstractArray{T}}) where {T} = zeroscalar(T)
zeroscalar(x::AbstractArray) = zeroscalar(eltype(x))

```

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [April 5, 2021, 6:19pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/65 "2021-04-05T18:19:15Z")

</div>

Split this topic on April 5, 2021, 6:19pm UTC

15 posts were split to a new topic: [Symbolic sqrt(2) and sin(pi)](https://discourse.julialang.org/t/symbolic-sqrt-2-and-sin-pi/58633)

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 4, 2021, 10:28pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/69 "2021-04-04T22:28:25Z")

</div>

```julia
julia> using Symbolics
julia> Term(sqrt,[2])
sqrt(2)

```

I don’t get why you’d work harder than that. You can build any expression to be lazily like that, or use `@register` on others.

> [@stevengj](#):
>
> This doesn’t work if `x` is an array — `zero(x) * x` will throw if `x` is not a square matrix. In the context of rings it is fine, of course, though it still won’t be equivalent to `x-x` for square matrices that may contain non-finite values.

I thought the reason for it was because the poster wanted “a zero in the type of `x`”:

> [@carette](#):
>
> I was warning that x−xx-x is tricky. It’s actually fairly safe to transform it to 0∗x0\*x . For one, that’s type-safe, in the sense that xx is still around, so you know the result’s type properly. (Polymorphic 00 is a major pain.)

In which case trying to get `0 * x` as they posted is an X-Y problem: just writing it as `zero(x)` is more direct and means what they actually meant using Base Julia functionality.

---

<div class="post-metadata">

**Author:** ![jzr](https://avatars.discourse-cdn.com/v4/letter/j/eb9ed0/32.png) [@jzr](https://discourse.julialang.org/u/jzr)\
**Post date:** [April 4, 2021, 11:16pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/70 "2021-04-04T23:16:50Z")

</div>

> [@ChrisRackauckas](#):
>
> I don’t get why you’d work harder than that. You can build any expression to be lazily like that, or use `@register` on others.

I think you’re right, I was mistaken. `Term(-, [x,x])` produces `x - x` which seems like what’s desired. Why doesn’t `x - x` produce the same result? For `@syms x::Any`, `x-x` errors, and for `@syms x::Num`, `x-x` returns `0::Float64`.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [April 5, 2021, 12:06am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/71 "2021-04-05T00:06:51Z")

</div>

> [@ChrisRackauckas](#):
>
> In which case trying to get `0 * x` as they posted is an X-Y problem: just writing it as `zero(x)` is more direct and means what they actually meant using Base Julia functionality.

`0 * x` is not the same as `zero(x)` if `x` is not finite (`NaN` or `±Inf`). That’s why `x - x` can be simplified to `0*x` (for an appropriate `0`) but not to `zero(x)` in general, and simplifying x-x is the [context](https://discourse.julialang.org/t/cas-best-practices/58092/61) here.

---

<div class="post-metadata">

**Author:** ![JeffreySarnoff](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jeffreysarnoff/32/1980_2.png) [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Post date:** [April 5, 2021, 10:34am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/74 "2021-04-05T10:34:03Z")

</div>

As the originator of this thread, it is my hope and request that we ( that means you Julians 🙂 ) return to being about “CAS Best Practices”. That there are different considerations held by CAS experts is helpful information. Let it be, so we move forward. Thank you all for participating.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 5, 2021, 10:39am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/75 "2021-04-05T10:39:07Z")

</div>

If you want it exact, don’t ask for the Float64 result?

```julia
julia> ex = Term(sqrt,[big(2)])
sqrt(2)

julia> eval(toexpr(ex))
1.414213562373095048801688724209698078569671875376948073176679737990732478462102

```

It just follows Julia’s semantics.

> [@stevengj](#):
>
> `0 * x` is not the same as `zero(x)` if `x` is not finite ( `NaN` or `±Inf` ). That’s why `x - x` can be simplified to `0*x` (for an appropriate `0` ) but not to `zero(x)` in general, and simplifying x−xx-x is the [context](https://discourse.julialang.org/t/cas-best-practices/58092/61) here.

I see. Yeah there’s a more interesting question to be answered here, which is what kind of number is the default symbol? The default is `@variables x::Real`, specifically not `::AbstractFloat`. Since these rules are all type-based, I think simplifying `x-x -> 0` on `::Real` is the semantics you’d expect, while `x-x -> 0 * x` should be what happens on `::AbstractFloat`. This keeps two algebras distinct: real numbers obey commutative and associative laws, they are a field, etc. while floating point numbers have a different set of rules.

The reason for keeping this difference is that the set of transformations that keep exact floating point semantics is already in LLVM, so I don’t think there’s too much to gain there. But it would be a good opt-in: “make simplifications only do what’s allowed on floating point numbers” instead of “make simplifications only do what’s allowed on real numbers”. The latter seems to be what people normally want on models, or at least out of a CAS, though there are times when exact floating point semantics need to be followed (in which case we have `@register` for someone to register a function and keep some calculation intact). What I want to do is add a call to [Herbie](https://herbie.uwplse.org/), so that in many cases we aren’t just applying rules that change floating point behavior but actually improve it, so even further down the path of “keep floating point separate”.

---

<div class="post-metadata">

**Author:** ![goretkin](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goretkin/32/167_2.png) [@goretkin](https://discourse.julialang.org/u/goretkin)\
**Post date:** [April 5, 2021, 7:33pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/76 "2021-04-05T19:33:43Z")

</div>

> [@ChrisRackauckas](#):
>
> The default is `@variables x::Real` , specifically not `::AbstractFloat` . Since these rules are all type-based, I think simplifying `x-x -> 0` on `::Real` is the semantics you’d expect, while `x-x -> 0 * x` should be what happens on `::AbstractFloat` .

This makes sense on its own, but I do not know how to square it with `AbstractFloat <: Real`. This is a question related to (or exactly about?) co/contra/invariance, and how the subtyping relation shows up, if at all, in figuring out what rules apply when. Is the following a bug? Or is it a necessary consequence of the rule defined on `Real`?

```julia
julia> y = SymbolicUtils.Sym{AbstractFloat}(:y)
y

julia> y - y
0

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 5, 2021, 7:40pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/77 "2021-04-05T19:40:50Z")

</div>

> [@goretkin](#):
>
> Is the following a bug? Or is it a necessary consequence of the rule defined on `Real` ?

I would say yes, it’s a missing specialization.

> [@goretkin](#):
>
> This makes sense on its own, but I do not know how to square it with `AbstractFloat <: Real` .

You have subtypes which override defaults of their supertype all of the time. `AbstractArray` for example has many array types which change underlying behaviors of broadcast, indexing, etc.

---

<div class="post-metadata">

**Author:** ![goretkin](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goretkin/32/167_2.png) [@goretkin](https://discourse.julialang.org/u/goretkin)\
**Post date:** [April 5, 2021, 9:46pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/78 "2021-04-05T21:46:19Z")

</div>

> [@ChrisRackauckas](#):
>
> You have subtypes which override defaults of their supertype all of the time. `AbstractArray` for example has many array types which change underlying behaviors of broadcast, indexing, etc.

But subtypes should not override the defaults (or “fallbacks”) of their supertype “willy nilly”. There should be some kind of consistency with the fallback. e.g.:

```julia
julia> typeof(1:4)
UnitRange{Int64}

julia> typeof(1:4) <: AbstractArray
true

julia> (1:4) .+ 5
6:9

julia> collect(1:4) .+ 5
4-element Vector{Int64}:
 6
 7
 8
 9

julia> collect((1:4) .+ 5) == collect(1:4) .+ 5 # consistency
true

```

If broadcasted addition on `UnitRange{Int64}` were defined in such a way that the “# consistency” property did not hold, that would be bad for generic programming!

Does that badness apply in this CAS setting? What consistency should exist across rules that are defined on types that have a subtype relationship among them? I am not sure. That’s what the question is.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 5, 2021, 9:56pm UTC](https://discourse.julialang.org/t/cas-best-practices/58092/79 "2021-04-05T21:56:57Z")

</div>

That’s a good question. I would assume that for number types you’d want to have those simplifications by default, and choose on a subset of `AbstractFloats` to disable them. I could see someone going the other way and having `@variables x` default to `SymbolicReal <: Real` and define those rules on `SymbolicReal`. It’s worth a discussion on the SymbolicUtils.jl repo. Open an issue.

---

<div class="post-metadata">

**Author:** ![sswatson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sswatson/32/5135_2.png) [@sswatson](https://discourse.julialang.org/u/sswatson)\
**Post date:** [April 6, 2021, 12:16am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/80 "2021-04-06T00:16:15Z")

</div>

> [@rfateman](#):
>
> It is not an identity that sqrt(x^2) =x. sqrt(9) is ±3.

There is no problem whatsoever with defining the square root function (`sqrt` or \sqrt{}) to be the principal square root. This is the standard approach, as described in the [Wikipedia article on square roots](https://en.wikipedia.org/wiki/Square_root) and as implemented in Julia, Python, Mathematica, Maple, Maxima, etc. Insisting that we define \sqrt{9} to be \pm3 is a pathway to confusion.

I would think that if you want to show that \sqrt{x^2} and x are unequal as expressions, you’d want to point out the example of x = -3. That way the argument does not depend on a non-standard interpretation of the square root symbol.

---

<div class="post-metadata">

**Author:** ![JeffreySarnoff](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jeffreysarnoff/32/1980_2.png) [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Post date:** [April 6, 2021, 12:27am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/81 "2021-04-06T00:27:30Z")

</div>

From a conventional programming language point-of-view, yes,

From a logico-symbolic point-of-view where  
`y = sqrt(x) means that x == y^2 `, `y = sqrt(9) means that 9 == y^2`  
asking Maple to `solve(y^2 = 9, y)` gives `(3, -3)`

A good deal of the back-and-forth is frisson between Numerics and Symbolics (both are math).

We are working toward a best (or a very good) case implementation approach that allows us to glean important benefits that each perspective brings, sometimes one better than the other and vice versa.

---

<div class="post-metadata">

**Author:** ![sswatson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sswatson/32/5135_2.png) [@sswatson](https://discourse.julialang.org/u/sswatson)\
**Post date:** [April 6, 2021, 4:20am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/82 "2021-04-06T04:20:36Z")

</div>

My point is exactly that the square root function is, in all common usage, distinct from the solve function. If we insist on conflating the two, how do we interpret an expression like \sqrt{2} + \sqrt{3}? It’s got to be a four-element set, I suppose. Here we’re already verging on chaos, wishing we had a concise way to say “the positive square root” when we want to.

I’m comfortable granting the point that there are many interesting loci of interplay between the demands of symbolic and numeric computation, some of which may require creative thinking about how common concepts need to be treated differently. But I remain skeptical about the specific proposal of abandoning the square root function (which is implied by an assertion like \sqrt{9} = \pm 3).

---

<div class="post-metadata">

**Author:** ![JeffreySarnoff](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jeffreysarnoff/32/1980_2.png) [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Post date:** [April 6, 2021, 4:51am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/83 "2021-04-06T04:51:36Z")

</div>

We are not going to abandon any of these `print(names(Base.Math))`. I understand how that impression may have been conveyed, inadvertently. What is being considered is how to respect the numerical math functions without abandoning the symbolic/theoretical math functions that may carry the same names in ordinary discourse (hence the tangling of notions). It is my guess that we will find a way to distinguish the desired context (symbolic, numeric, some of each) for any occurrence of the named function – there are multiple avenues for this. The simplest way to distinguish the Julia numerical square\_root from the more constructively mathematical symbolic square\_root would be use the name from Base.Math as the numerical function and to use its titlecase spelling as the symbolic function `sqrt(2) !== Sqrt(2)`. That leaves out a mixed version (maybe `Sqrt` is mixed, and purely symbolic is a context, as with `let`, and `do` blocks).

see next post by @DNF for one of the more familiar approaches within Julia

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [April 6, 2021, 7:44am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/84 "2021-04-06T07:44:44Z")

</div>

> [@JeffreySarnoff](#):
>
> use the name from Base.Math as the numerical function and to use its titlecase spelling as the symbolic function `sqrt(2) !== Sqrt(2)`

Doesn’t that run a bit counter to the general Julia philosophy? Could one instead leverage multiple dispatch with `sqrt(2)` vs `sqrt(Symbolic(2))` (just making up names here)?

Disclaimer: I know nothing about CASs in general, and hardly anything about Symbolics.jl, so this is a question rather than a suggestion.

---

<div class="post-metadata">

**Author:** ![JeffreySarnoff](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jeffreysarnoff/32/1980_2.png) [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Post date:** [April 6, 2021, 8:19am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/85 "2021-04-06T08:19:48Z")

</div>

of course. I just wanted to clarify that these two things (the numerical evaluation and the symbolic handling) are not required to be one thing (although synthesis is good if it work).

---

<div class="post-metadata">

**Author:** ![Fredrik\_Johansson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fredrik_johansson/32/23679_2.png) [@Fredrik\_Johansson](https://discourse.julialang.org/u/Fredrik_Johansson)\
**Post date:** [April 6, 2021, 8:43am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/86 "2021-04-06T08:43:05Z")

</div>

The principal square root function is a logically, mathematically, symbolically and numerically perfectly fine function. Using the names sqrt/Sqrt to denote anything else is a mistake. Not because other notions of square root aren’t useful too, but because this convention is already firmly established. If you want a function that returns a set, call it something else (all\_square\_roots(x) or whatever).

The more general issue here is that computer algebra systems historically are very confused about what it means to solve something. Does solve() return some particular solution, the set of solutions, a set containing some of the solutions, a superset of the solutions (possibly containing spurious solutions), real solutions, complex solutions, or something random depending on context or the phase of the moon? I’ve seen all combinations.

This confusion essentially goes away if you build symbolic sets with well-typed quantifiers into the system. If you want to solve x^2 = y, write down S = \{x : x \in \mathbb{C} \land x^2 = y \} or S = \{x : x \in \mathbb{R} \land x \ge 0 \land x^2 = y \} or something else, depending on what you mean. Now various notions of solving can be reduced to computations on the set S: checking its cardinality, extracting an arbitrary element, iterating over the elements, computing some superset, taking the intersection with another set, etc. To keep things convenient for users, you may want unique\_solution(), any\_solution() and all\_solutions() functions or something similar.

There is actually nothing inherently wrong with simplifying sqrt(x^2) → x or solving x^2 = y → x = sqrt(x) _as purely formal operations on formal symbolic expressions_. You just have to remember that these are formal operations which don’t commute with evaluation for x \in \mathbb{R}, say. When we do mathematics by hand, we often ignore conditions and check that the result makes sense in the end. The traditional CAS design approach mimics human operation but leaves the final checking to the user. This is clearly very useful in many cases, but it makes it seriously difficult to write correct complex programs that compose operations (where the user can’t realistically check each step). I’m convinced that it’s possible to do better; there’s no need to throw out the classical heuristic symbolic manipulation completely, but it should be a complement to and not a substitute for correct domain-aware rewriting. There are some parallels with floating-point versus interval arithmetic here (including the fact that floating-point solutions much like heuristic symbolic solutions sometimes can be checked a posteriori!).

I also completely agree with Oscar Benjamin that automatic evaluation is a mistake. Simple, unevaluated expressions by default are a win-win: easier to implement, easier and more versatile to use (the user can decide precisely which transformations to apply), and more efficient (symbolic expressions are a terrible data structure for algebraic computation). Instead of performing clever and invariably buggy automatic “canonicalization” on expressions, a far better solution is to convert to specialized data structures internally for computations (e.g. polynomial types for expanding or factoring).

---

<div class="post-metadata">

**Author:** ![JeffreySarnoff](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jeffreysarnoff/32/1980_2.png) [@JeffreySarnoff](https://discourse.julialang.org/u/JeffreySarnoff)\
**Post date:** [April 6, 2021, 8:45am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/87 "2021-04-06T08:45:37Z")

</div>

Thank you @Fredrik_Johansson. Your comments are clear and on point.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 6, 2021, 10:15am UTC](https://discourse.julialang.org/t/cas-best-practices/58092/88 "2021-04-06T10:15:08Z")

</div>

Yes, that’s why in the tracer it’s planned to not do automatic evaluation, as described in [Symbolic sqrt(2) and sin(pi) - #21 by ChrisRackauckas](https://discourse.julialang.org/t/symbolic-sqrt-2-and-sin-pi/58633/21), with dispatches for alternative constructor-based simplifications all being type based.

> [@Fredrik\_Johansson](#):
>
> To keep things convenient for users, you may want unique\_solution(), any\_solution() and all\_solutions() functions or something similar.

That’s a good idea.

> [@Fredrik\_Johansson](#):
>
> Simple, unevaluated expressions by default are a win-win: easier to implement, easier and more versatile to use (the user can decide precisely which transformations to apply), and more efficient (symbolic expressions are a terrible data structure for algebraic computation). Instead of performing clever and invariably buggy automatic “canonicalization” on expressions, a far better solution is to convert to specialized data structures internally for computations (e.g. polynomial types for expanding or factoring).

That is what we do. Sounds like Symbolics.jl is on the right path then, thanks.

[Previous page](https://discourse.julialang.org/t/cas-best-practices/58092.md?page=3)

[Next page](https://discourse.julialang.org/t/cas-best-practices/58092.md?page=5)