Add `lambdify` or `build_function` for converting `GiacExpr` into native Julia callables

[s-celles]

Ref: https://discourse.julialang.org/t/ann-giac-jl-julia-interface-to-the-giac-computer-algebra-system/136681/6

Summary

It would be useful to have a built-in helper that turns a GiacExpr into a Julia function — both as an ergonomic shortcut over substitute + to_julia, and (optionally) as a way to compile the expression into native Julia code that no longer round-trips through Giac on each call.

(AI generated content below)

This is directly analogous to two well-established APIs in the Julia/Python symbolic ecosystem:

• SymPy's lambdify — converts a SymPy expression into a fast numeric function. Available in Julia via SymPy.jl's lambdify, which wraps SymPy's implementation. Users can write f = lambdify(expr, (x, y)) and get back a callable that evaluates the expression numerically.
• Symbolics.jl's build_function — generates native Julia code from a symbolic expression, optionally compiling in-place / array-output / GPU-targeted variants. This is the SciML-native equivalent and powers a lot of ModelingToolkit.jl's runtime performance.

Giac.jl currently has no equivalent. Users have to manually compose substitute + to_julia, which works but is undiscoverable for newcomers coming from either ecosystem.

Motivation

The current idiom for plotting or numerical evaluation is:

@giac_var x y a b
expr = sin(√(a * x^2 + b * y^2)) * cos(x/2)

num_expr = substitute(expr, Dict(a => 1.0, b => 1.0))
f(_x, _y) = to_julia(substitute(num_expr, Dict(x => _x, y => _y)))

# or with v0.12.0
num_expr = expr(a => 1.0, b => 1.0)

f(_x, _y) = to_julia(num_expr(x => _x, y => _y))

This works, and the v0.12 improvements — call-syntax, simultaneous semantics, to_julia auto-evalf, direct giac_subst binding — already make it cleaner and faster than v0.11. But:

1. Users coming from SymPy / SymPy.jl will look for lambdify and not find it. Users coming from Symbolics.jl will look for build_function and not find it. A single named entry point closes that gap.
2. Every call still crosses the Julia ↔ C++ boundary, which is fine for a 100×100 surface plot (~10k calls) but becomes a bottleneck for interactive Pluto sliders on larger grids, parameter sweeps, or numerical optimization where f is called millions of times. SymPy's lambdify and Symbolics' build_function both solve exactly this problem by emitting native code.

A lambdify helper would document the recommended pattern and give a clean upgrade path to a compiled version when performance matters.

Proposed API

Tier 1 — convenience wrapper (always available, stays inside Giac)

Roughly equivalent to SymPy's default lambdify behavior, but staying within Giac for evaluation:

"""
    lambdify(expr::GiacExpr, vars::GiacExpr...) -> Function

Return a Julia function `f(vals...)` that evaluates `expr` by
simultaneously substituting `vars` with `vals` and converting
the result via `to_julia`.

Mirrors `SymPy.lambdify` and `Symbolics.build_function` in spirit:
turn a symbolic expression into a numeric callable.
"""
lambdify(expr, vars...) =
    (vals...) -> to_julia(expr(Pair.(vars, vals)...))

This is essentially documentation + a single exported symbol. Cheap to implement, gives users the API they expect.

Tier 2 — compiled version (opt-in, leverages v0.12 introspection)

This is what makes SymPy's lambdify and Symbolics' build_function actually fast: emit native code instead of calling back into the CAS for each evaluation.

The new v0.12 introspection helpers (iscall, operation, arguments, free_symbols, unwrap_const, plus the GiacTermInterfaceExt weak-dep extension) make this directly feasible — they expose exactly the AST-walking interface needed:

lambdify(expr, vars...; compile::Bool = false)

When compile=true, the implementation walks the AST, maps Giac heads to Julia equivalents (sin → sin, + → +, …), and @evals an anonymous function. Subsequent calls run as plain Julia code — no Giac round-trip, LLVM-vectorizable.

Sketch of the AST walk:

function _giac_to_expr(g, lookup)
    if !iscall(g)
        c = unwrap_const(g)
        return c === nothing ? Symbol(string(g)) : c
    end
    op   = Symbol(string(operation(g)))
    args = _giac_to_expr.(arguments(g), Ref(lookup))
    jop  = get(lookup, op, op)
    return Expr(:call, jop, args...)
end

A GIAC_TO_JULIA lookup table handles cases where Giac and Julia names diverge. SymPy maintains a similar table (MODULES in sympy.utilities.lambdify) and Symbolics.jl handles this through its toexpr machinery — both are good references for what coverage looks like in practice.

Optional Tier 3 — bridge to Symbolics.jl

Rather than reimplementing the full toolchain, a third backend could simply delegate to Symbolics.build_function via the existing to_symbolics extension:

lambdify(expr, vars...; backend = :giac)
# backend ∈ (:giac, :compile, :symbolics)

This gives Giac.jl users access to everything Symbolics.jl already does well — in-place mutation, array outputs, sparsity, cse, GPU codegen via @target_func, etc. — without Giac.jl needing to maintain any of that infrastructure. It's the same pragmatic choice SymPy.jl makes by deferring to upstream SymPy's lambdify.

Open questions

1. Scalar-only or vector-aware? Both SymPy's lambdify and Symbolics' build_function support array outputs and in-place forms (build_function(expr, args...; target=Symbolics.JuliaTarget()) returns both an out-of-place and in-place version). Worth scoping for v1, or punt to a follow-up?
2. Compiled-mode coverage. What's the right minimum set of Giac heads to support out of the box, and what should the fallback be when an unknown head is encountered (error vs. fall back to Tier 1 evaluation of that subexpression)? SymPy errors loudly; Symbolics tends to pass through.
3. World-age. Tier 2 uses @eval, which means an immediate call in the same global scope can hit the world-age boundary. SymPy.jl works around this with Base.invokelatest; worth documenting either way. For plotting use cases (the most common one) it's a non-issue because Plots calls f from inside its own functions.
4. Constants. How does unwrap_const currently expose π, ℯ, ∞? They need to map to Base.pi, Base.ℯ, Inf in the emitted code (SymPy's lambdify does this via its MODULES mapping).
5. Naming. lambdify is the SymPy name and is widely recognized in both Python and Julia (via SymPy.jl). build_function is the Symbolics.jl name. Going with lambdify probably maximizes discoverability — users coming from either ecosystem will recognize the concept, and the name is shorter and unambiguously refers to numeric evaluation. Symbolics users typing build_function could be helped by a docstring @doc alias or just a mention in the docs.

Proposed scope for a first PR

Just Tier 1, exported as lambdify, plus a docs page showing the three patterns (manual substitute, lambdify, lambdify(...; compile=true)) side by side, with a comparison table to SymPy's lambdify and Symbolics' build_function so users from either ecosystem can map their mental model onto Giac.jl quickly. Tier 2 and Tier 3 can land as follow-ups once the API surface is agreed upon.

---

Trade-offs of each tier

Tier 1 — convenience wrapper (stays inside Giac)

The deal: lambdify becomes a one-liner that wraps the existing substitute + to_julia pattern. Each call still crosses into Giac/C++.

Strengths

• Zero new dependencies. Works the moment a user types using Giac. No Symbolics.jl, no codegen machinery.
• Trivial to implement and maintain. ~5 lines of code; nothing to break.
• Zero divergence risk. Whatever Giac computes is what the user gets — no parallel evaluator that could drift from Giac's semantics.
• Handles every Giac expression. If substitute and to_julia work on it, lambdify works on it. No missing operators, no AST coverage gaps.
• Predictable for newcomers. The mental model is "it's just substitute with a nicer name."

Weaknesses

• Slow per-call. Every evaluation pays the Julia ↔ C++ FFI cost: marshalling Julia values into Giac generic objects, running Giac's substitution and simplification, marshalling back. Order-of-magnitude estimate: ~1–10 µs/call rather than ~10–100 ns for native Julia. Fine for 10k-point surface plots, painful for 10M-point optimization loops.
• Not autodiff-friendly. ForwardDiff.gradient(f, x) won't work because Dual numbers can't traverse the Giac boundary. Same for Zygote, Enzyme.
• No vectorization. LLVM can't see through the FFI call, so SIMD and loop unrolling are off the table.
• Allocations. Each call typically allocates at least one GiacExpr for the substitution result.

Best for: plots, exploratory analysis, teaching, anywhere call count is modest and the value is in the convenience.

---

Tier 2 — compile to native Julia (in-house codegen)

The deal: Walk the GiacExpr AST, emit a Julia Expr, @eval it. Subsequent calls are pure Julia, JIT-compiled by LLVM, no Giac involvement.

Strengths

• Native speed. Once compiled, performance is identical to a hand-written Julia function. SIMD-vectorizable, inlinable, unboxed.
• Autodiff-compatible. The emitted function is plain Julia operating on whatever number type comes in, so ForwardDiff and friends work transparently — provided the operator table maps to generic functions like Base.sin rather than Float64-specialized variants.
• Zero runtime dependency on Giac. After compilation, the function survives even if the Giac library is unloaded.
• No new package dependency. Lives entirely inside Giac.jl, using only the v0.12 introspection helpers (iscall, operation, arguments, unwrap_const).

Weaknesses

• Non-trivial implementation surface. The AST walker plus the Giac→Julia operator mapping table must cover everything users throw at it. Edge cases include: rationals (3/4 should stay rational, not become 0.75 prematurely), arbitrary-precision integers, π/ℯ/∞, branch cuts (does Giac's log agree with Julia's on negative reals?), piecewise functions, unevaluated integrals/sums/limits, user-defined functions, Ei/Si/erf and other special functions that Giac knows but Base doesn't.
• Maintenance burden grows over time. Each new Giac function someone wants to use forces a table update.
• Divergence risk. If the operator mapping is wrong (e.g. Giac's log is natural log but the table maps to log10), the compiled function silently produces wrong numbers. Unit tests must compare against Tier 1 evaluation on a sample grid.
• World-age friction. @eval creates a function in a newer world than the caller's. Calling it immediately from the same global scope hits MethodError: no method matching … (method too new). Mitigations: Base.invokelatest, or document that compiled functions should be passed to consumers (like Plots.surface) rather than called directly at top level.
• Compilation latency. The first call to a freshly-compiled function pays Julia's TTFX cost. For interactive use with frequently-changing expressions (e.g. Pluto sliders that re-trigger compilation on each parameter change), this can be more annoying than Tier 1's per-call FFI.
• Allocation if not careful. Naive codegen can produce intermediate temporaries; well-tuned codegen avoids them. Symbolics.jl spent years on this.

Best for: hot loops where the same expression is evaluated millions of times with fixed structure, autodiff-driven workflows, performance-critical applications where the user is willing to pay the upfront cost.

---

Tier 3 — delegate to Symbolics.jl build_function

The deal: Convert GiacExpr → Num via the existing to_symbolics extension, then hand off to Symbolics.build_function, which already knows how to do everything Tier 2 promises.

Strengths

• Industrial-strength codegen, for free. Symbolics.jl has years of tuning behind build_function: out-of-place + in-place pairs, common subexpression elimination, sparsity detection, expression=true for inspectable code, alternative targets (CTarget, StanTarget, GPU via the SciML stack).
• Autodiff-native. The resulting function is plain Julia and integrates cleanly with the entire SciML ecosystem (DifferentialEquations.jl, Optimization.jl, ModelingToolkit.jl).
• Tiny implementation surface in Giac.jl. ~10 lines in a new weak-dep extension. The hard work happens upstream.
• Operator coverage tracks Symbolics.jl. When Symbolics.jl learns a new function, Giac.jl gets it for free — no table to update, no PR to chase.
• Vector-of-expressions for free. lambdify([expr1, expr2, expr3], (x, y); backend=:symbolics) returns a tuple (f_oop, f_ip!) — exactly what RHS-of-ODE workflows need.

Weaknesses

• Heavy optional dependency. Symbolics.jl plus its transitive closure (SymbolicUtils, TermInterface, Setfield, RuntimeGeneratedFunctions, …) adds substantial precompile time. Educational users on slower machines may notice. Mitigated by the weak-dep extension model — the cost is paid only when using Symbolics is loaded.
• Round-trip lossiness. to_symbolics ∘ to_giac is not guaranteed to be the identity. Giac may simplify in one direction, Symbolics in another. Most of the time this is benign (both agree numerically), but it can surface in:
  • Branch choices. √(x²) vs |x| vs x.
  • Domain assumptions. Giac may assume real positive variables in some commands; Symbolics treats variables as complex by default.
  • Special functions. Functions Giac knows but Symbolics doesn't (or vice versa) hit the boundary as opaque atoms.
• Audit needed. Before this tier can be advertised as "just works," to_symbolics needs systematic testing on the full range of Giac heads — trig, exp/log, roots, piecewise (when), unevaluated integrals/sums/limits, Heaviside, Dirac. Bugs here show up as wrong numerical results downstream.
• Two indirection layers between user and result. Debugging a bad output means asking: did Giac compute the wrong thing? Did to_symbolics translate it wrongly? Did build_function codegen it wrongly? More moving parts means more places to look.
• API stability coupling. If Symbolics.build_function ever changes its signature or default options, Giac.jl's :symbolics backend has to track. Manageable, but a real cost.

Best for: users already working in the SciML ecosystem, anyone who needs in-place forms / array outputs / autodiff / sparsity, advanced workflows where the dependency cost is dwarfed by the capabilities gained.

---

Summary table

Concern	Tier 1 (Giac)	Tier 2 (compile)	Tier 3 (Symbolics)
New dependency	None	None	Symbolics.jl (weak)
LOC in Giac.jl	~5	~50–150	~10
Per-call cost	Giac FFI (~µs)	Native Julia (~ns)	Native Julia (~ns)
First-call latency	~zero	Compile time	Compile time + load
Autodiff	✗	✓ (if table is generic)	✓
Array outputs	Manual	DIY	✓ free
In-place forms	✗	DIY	✓ free
CSE	✗	DIY	✓ free
GPU/alt targets	✗	DIY	✓ free
Operator coverage	All of Giac	Whatever the table covers	Whatever to_symbolics covers
Divergence risk	None	Operator-table bugs	Round-trip semantics
Maintenance	Trivial	Ongoing (table)	Trivial (delegated)
Best for	Plots, teaching, exploration	Hot loops, embedded use, no-deps perf	SciML workflows, autodiff, vector outputs

---

How they compose

These tiers aren't mutually exclusive — they layer cleanly:

• Tier 1 is the floor. It always works, even with no extras loaded. Default backend = :giac.
• Tier 2 is for users who want native speed without pulling Symbolics.jl. Useful for educational settings, embedded contexts, or anyone wary of dependency weight.
• Tier 3 is for users already in the SciML world or who need features (in-place, CSE, autodiff, sparsity) that would be wasteful to reimplement.

A reasonable shipping order is Tier 1 first (lands the API), Tier 3 second (cheap to add, big payoff for power users), Tier 2 third or never (only if there's demonstrated demand for native speed without Symbolics.jl). Tier 2 is the most work for the narrowest audience — it might honestly never be worth building if Tier 3 covers the performance use cases adequately.

Reference:

• build_function https://symbolics.juliasymbolics.org/dev/manual/build_function/
• lambdify https://docs.sympy.org/latest/modules/utilities/lambdify.html

[s-celles]

Tier 1 & 3 implemented
Tier 2 deferred

[s-celles]

@tduretz Good news. I added a build_function that does exactly what you asked about . It's the direct counterpart to SymPy's lambdify / julia_code, and the name follows the Symbolics.jl / ModelingToolkit.jl convention.

Two backends are available via the backend keyword:

• backend = :giac (default) - a thin wrapper that stays inside Giac. The returned function does substitute + to_julia on each call. No extra dependency, covers every Giac expression, great for plotting and interactive exploration.
• backend = :symbolics converts the expression via to_symbolics and delegates to Symbolics.build_function. You get the full SciML machinery: native Julia code (~ns/call instead of ~µs through the FFI), out-of-place and in-place forms, autodiff support (ForwardDiff, Zygote, Enzyme), CSE (Common Subexpression Elimination), array outputs, etc. It's wired through a weak-dep extension, so the cost is only paid when using Symbolics is loaded.

Example:

using Giac
@giac_var x y
expr = sin(x^2 + y^2)

# Giac backend (default) — always available
f = build_function(expr, x, y)
f(1.0, 2.0)

# Symbolics backend — native Julia, autodiff-compatible
using Symbolics
f_fast = build_function(expr, x, y; backend = :symbolics)
f_fast(1.0, 2.0)

For inspectable Julia code (the julia_code() analogue), the :symbolics backend also gives you access to expression = Val{true} from Symbolics.build_function.

Design notes, backend trade-offs, and discussion of a possible third tier (tier 2 ie native codegen with no Symbolics.jl dependency) are in this issue but have been deferred