Questions about Simplification and Performance in Symbolics.jl

General Usage Performance
1

Dear all,

I have two questions regarding the use of Symbolics.jl, particularly around simplification and performance behavior.

1. Mathematical Identity Simplification

It seems that Symbolics.jl does not simplify certain mathematical identities automatically. For example:

using Symbolics
using BenchmarkTools

@variables x y

f(x, y) = x^2 + y^2 + 2x*y
f_expr = f(x, y)
f_simplified = Symbolics.simplify(f_expr)

This returns:

f_expr = f(x, y) 

x^2 + 2x*y + y^2

f_simplified = Symbolics.simplify(f_expr) 

x^2 + 2x*y + y^2

In this case, it does not simplify the expression to the expected identity form:(x + y)^2

I understand that SymbolicUtils.jl and MetaTheory.jl are used as backends. However, I’m wondering:
Do we need to manually define rules for standard mathematical identities, such as the binomial expansion, or is there a more built-in approach for this?

2. Performance & Allocation with build_function

Here’s another behavior I found interesting and a bit confusing related to performance after simplification.

Consider the following:

using Symbolics
using BenchmarkTools

@variables x y

f(x, y) = x + y^2 + x + y^2

@show f_expr = f(x, y)
@show f_simplified = Symbolics.simplify(f_expr)

f_built = build_function(f_simplified, x, y)
f_optimized = eval(f_built)

This produces:

f_expr = f(x, y) = 2x + 2(y^2)
2x + 2(y^2)

f_simplified = Symbolics.simplify(f_expr) = 2(x + y^2)
2(x + y^2)

So far, so good , the expression is simplified correctly.
However, performance differs significantly:

x_, y_ = 1.0, 2.0

@btime f($x_, $y_)  # original function
# 1.400 ns (0 allocations: 0 bytes)

@btime f_optimized($x_, $y_)  # built from Symbolics
# 18.637 ns (3 allocations: 48 bytes)

Looking at the LLVM output explains some of this:

Original function:

@code_llvm f(x_, y_)

Generates inline fmul instructions, like:

; Function Signature: f(Float64, Float64)
;  @ c:\Users\lab30501\Desktop\test\MWE.jl:6 within `f`
; Function Attrs: uwtable
define double @julia_f_25777(double %"x::Float64", double %"y::Float64") #0 {
top:
; ┌ @ intfuncs.jl:370 within `literal_pow`
; │┌ @ float.jl:493 within `*`
    %0 = fmul double %"y::Float64", %"y::Float64"
; └└
; ┌ @ operators.jl:596 within `+` @ float.jl:491
   %1 = fadd double %0, %"x::Float64"
   %2 = fadd double %1, %"x::Float64"
; │ @ operators.jl:596 within `+`
; │┌ @ operators.jl:553 within `afoldl`
; ││┌ @ float.jl:491 within `+`
     %3 = fadd double %0, %2
     ret double %3
; └└└
}

Optimized (Symbolics-generated) function:

@code_llvm f_optimized(x_, y_)

calls @j_pow:

; Function Signature: var"#13"(Float64, Float64)
;  @ none within `#13`
; Function Attrs: uwtable
define double @"julia_#13_25782"(double %"x::Float64", double %"y::Float64") #0 {
top:
;  @ none within `#13` @ C:\Users\lab30501\.julia\packages\Symbolics\BefTf\src\build_function.jl:137
; ┌ @ C:\Users\lab30501\.julia\packages\SymbolicUtils\aooYZ\src\code.jl:411 within `macro expansion`
; │┌ @ math.jl:1199 within `^`
    %0 = call double @j_pow_body_25786(double %"y::Float64", i64 signext 2)
; │└
; │┌ @ float.jl:491 within `+`
    %1 = fadd double %0, %"x::Float64"
; │└
; │┌ @ promotion.jl:430 within `*` @ float.jl:493
    %2 = fmul double %1, 2.000000e+00
    ret double %2
; └└
}

So it seems that instead of using y*y, the symbolic version uses a general-purpose power function, even for small integer powers like 2. This leads to additional function calls and allocations.

I’d really appreciate any insights or suggestions to understand these behaviors better. Thank you!

1 Like