JuliaSymbolics/Symbolics.jl

Simplified derivative expansion gives float coefficients

hersle opened this issue · 7 comments

The simplified derivative expansion

using Symbolics
@variables t a(t)
D(expr) = expand_derivatives(Differential(t)(expr))
simplify_fractions(D(D(a)/a) + D(D(a)/a))

outputs

(2.0a(t)*Differential(t)(Differential(t)(a(t))) - 2.0(Differential(t)(a(t))^2)) / (a(t)^2)

Why are the coefficients 2.0 instead of 2? I would like them to be exact integers/rationals instead of floats. Thanks!

It's a bug, we should find how what causes it.

Here is the most minimal example I am able to reduce it to:

using Symbolics
@variables a b
simplify(a/b + 2*a^2/b)

This outputs (a + 2.0(a^2)) / b, but I want (a + 2(a^2)) / b. Replacing simplify with simplify_fractions changes nothing.

With expand(a/b + 2*a^2/b) instead, I get a / b + (2(a^2)) / b, so perhaps the problem is related to how simplify puts everything in one fraction?

Possibly, and then it evaluates the fraction with / which gives a float back for two ints.

I've done some digging. It seems that even

using SymbolicUtils
@syms x
pf = PolyForm(x) # polynomial 1*x^1
div(pf, pf) # x / x

outputs 1.0 instead of 1.

Digging into in SymbolicUtils.jl's polyform.jl, it seems that it ultimately relies on DynamicPolynomials.jl to do the division.

There is an identical open/unresolved issue in that library. But it seems to me that this SymbolicUtils.jl commit was made to provide a workaround?

What is the way to go about this? I suppose the ideal way would be to fix it upstream?

All examples in this issue would be fixed by JuliaAlgebra/MultivariatePolynomials.jl#294.

This was fixed by JuliaAlgebra/MultivariatePolynomials.jl/pull/296 😃

using Symbolics
@variables t a(t)
D(expr) = expand_derivatives(Differential(t)(expr))
simplify_fractions(D(D(a)/a) + D(D(a)/a))

now gives rational coefficients:

((2//1)*a(t)*Differential(t)(Differential(t)(a(t))) - (2//1)*(Differential(t)(a(t))^2)) / (a(t)^2)

Similarly,

using Symbolics
@variables a b
simplify(a/b + 2*a^2/b)

gives

(a + (2//1)*(a^2)) / b