` exploring – which is how I do most of my package capabilities exploring.
* * *
A related (and perhaps more uplifting) comment, I’m coming back to Julia now to convert some educational material from Matlab and the performance between Julia and Matlab symbolic calculations is ridiculous. Obviously, Matlab supports a broader range of functionality (e.g. multivariate integration), but it appears Julia has a good foundation. Here’s a direct (manual) translation of Matlab code into Julia. Since the ultimate goal is educational, I’m trying to avoid too many performance tricks that might obfuscate the code for non-CS undergrads.
```julia
function LagrangeBasis( degree, basis_idx, variate, domain )
ref_domain = PolynomialReferenceDomain( Lagrange )
variate = AffineMapping1D( domain, ref_domain, variate )
node = LinRange( -1.0, 1.0, degree + 1 )
basis_val = 1.0
for n = 1 : degree + 1
if n != basis_idx
basis_val = basis_val * ( (variate - node[n] ) / ( node[basis_idx] - node[n] ) )
end
end
return basis_val
end
function PolynomialReferenceDomain( basis::Basis )
if basis == Bernstein
return [0.0, 1.0]
elseif basis == Chebyshev
return [-1.0, 1.0]
elseif basis == Lagrange
return [-1.0, 1.0]
elseif basis == LagrangeCheby
return [-1.0, 1.0]
elseif basis == Legendre
return [-1.0, 1.0]
elseif basis == Monomial
return [0.0, 1.0]
end
end
function AffineMapping1D( domain, target_domain, x )
x_target = x .- domain[1]
x_target *= ( ( target_domain[2] - target_domain[1] ) / ( domain[2] - domain[1] ) )
x_target += target_domain[1]
end
@enum Basis begin
Bernstein
Chebyshev
Lagrange
LagrangeCheby
Legendre
Monomial
end
```
Benchmarking:
```julia
julia> degree = 8
julia> basis_idx = 7
julia> domain = [0, 1]
julia> BenchmarkTools.@btime( LagrangeBasis( $degree, $basis_idx, $x, $domain ) )
34.200 μs (440 allocations: 21.95 KiB)
91.02222222222221x*(2.0x - 0.25)*(2.0x - 0.5)*(2.0x - 0.75)*(2.0x - 1)*(2.0x - 1.25)*(2.0x - 1.75)*(2.0x - 2.0)
```
Whereas Matlab’s `timeit` function on its equivalent code (using Symbolic Toolbox calculations) results in `0.0196 sec` execution time – meaning Julia’s `Symbolics.jl` is roughly 500x-600x faster (in this instance) without the user really worrying about performance.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [August 7, 2023, 4:16pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/4 "2023-08-07T16:16:56Z")
Yeah I think that’s a good way to put it. Symbolics has a good foundation but needs more features to fill it out. That is happening over time, but integration is one piece that’s notably missing.
---
### Author:  [@Ashwani\_Rathee](https://discourse.julialang.org/u/Ashwani_Rathee)
#### Post date: [August 7, 2023, 6:18pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/5 "2023-08-07T18:18:54Z")
Does this also affect MethodOfLines.jl directly which tells it doesn’t support IntegralEquations for discretization? given Symbolics.jl is used for defining
---
### Author:  [@MilesCranmer](https://discourse.julialang.org/u/MilesCranmer)
#### Post date: [August 7, 2023, 7:40pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/6 "2023-08-07T19:40:50Z")
If needed, you can always fall back to SymbolicRegression.jl and **learn** the symbolic integral with a genetic algorithm, by evaluating the _derivatives_ of candidate expressions and seeing if they match the function. (This is actually an even more powerful strategy than available in Matlab/Mathematica/SymPy!)
Then you can just convert the recovered expression to Symbolics.jl afterwards with `node_to_symbolic`!
Here’s a demo of integrating the function f(x) = x^2 + 2x + 1 + \sin(x) by numerically evaluating its derivatives inside a genetic algorithm:
[https://global.discourse-cdn.com/julialang/original/3X/c/4/c404b18583793395d75852f6b8b2ce79f4ac438c.mp4](https://global.discourse-cdn.com/julialang/original/3X/c/4/c404b18583793395d75852f6b8b2ce79f4ac438c.mp4)
Code:
```julia
using MLJ # for fit/predict
using SymbolicRegression # for SRRegressor
using Zygote # For `enable_autodiff=true`
using SymbolicUtils
f(x) = x^2 + 2x + 1 + sin(x)
X = reshape(0.0:0.01:10.0, :, 1)
y = f.(X[:, 1])
function derivative_loss(tree, dataset::Dataset{T,L}, options, idx) where {T,L}
# Column-major:
X = idx === nothing ? dataset.X : view(dataset.X, :, idx)
∂y = idx === nothing ? dataset.y : view(dataset.y, idx)
ŷ, ∂ŷ, completed = eval_grad_tree_array(tree, X, options; variable=true)
!completed && return L(Inf)
mse = sum(i -> (∂ŷ[i] - ∂y[i])^2, eachindex(∂y)) / length(∂y)
return mse
end
model = SRRegressor(
binary_operators=[+, -, *, /],
unary_operators=[sin, cos, exp],
loss_function=derivative_loss,
enable_autodiff=true,
batching=true,
batch_size=25,
niterations=100,
)
mach = machine(model, X, y)
fit!(mach)
r = report(mach)
eq = r.equations[r.best_idx]
symbolic_eq = node_to_symbolic(eq, model)
```
---
### Author:  [@GregVernon](https://discourse.julialang.org/u/GregVernon)
#### Post date: [August 7, 2023, 7:58pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/7 "2023-08-07T19:58:48Z")
While certainly (very) interesting and (very) cool, it’s a little heavy for the purposes of teaching undergraduates L^2 inner products:
> Let’s compute the inner product between x^0 and x^1 on the domain [-1,1]. First we will train a symbolic regression model to integrate the integrand with a genetic algorithm – an NP-hard problem…
Meanwhile the students be like:
> Why all this machinery when \int x \ \mathrm{dx} is obviously just \frac{x^2}{2} ?
It’s just much cleaner right now in my Matlab-based materials:
```julia
int( x^0 * x^1, [-1, 1] ) % --> 0
```
---
### Author:  [@MilesCranmer](https://discourse.julialang.org/u/MilesCranmer)
#### Post date: [August 7, 2023, 8:20pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/8 "2023-08-07T20:20:58Z")
Ah, in that case I misread 😂
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [August 7, 2023, 10:35pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/9 "2023-08-07T22:35:58Z")
> [@Ashwani\_Rathee](#):
>
> Does this also affect MethodOfLines.jl directly which tells it doesn’t support IntegralEquations for discretization? given Symbolics.jl is used for defining
No, the ability to represent integrals is separate from the ability to symbolically solve them. Symbolics has the ability to represent integrals right now, which can be used for example to represent integro-differential equations for NeuralPDE and MOL, but it cannot symbolically solve integrals.
> [@MilesCranmer](#):
>
> If needed, you can always fall back to SymbolicRegression.jl and **learn** the symbolic integral with a genetic algorithm, by evaluating the _derivatives_ of candidate expressions and seeing if they match the function. (This is actually an even more powerful strategy than available in Matlab/Mathematica/SymPy!)
You may want to check out SymbolicNumericIntegration.jl which could in theory have a SymbolicRegression.jl backend to how it’s doing the symbolic regression part. It’s a bit smarter than just a direct symbolic regression though, I think you’d enjoy it.
> **[Symbolic-Numeric Integration of Univariate Expressions based on Sparse...](https://arxiv.org/abs/2201.12468)**
>
> Most computer algebra systems (CAS) support symbolic integration as core functionality. The majority of the integration packages use a combination of heuristic algebraic and rule-based (integration table) methods. In this paper, we present a hybrid...
But indeed it’s the kind of algorithm that should exist in a suite of symbolic integration algorithms, and this probably isn’t the the kind of thing you’d tell people to use first. RUBI type methods should come first.
---
### Author:  [@GregVernon](https://discourse.julialang.org/u/GregVernon)
#### Post date: [September 20, 2024, 6:48am UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/10 "2024-09-20T06:48:46Z")
Well, I’m back at it again – needing to integrate Symbolic expressions as part of some educational collateral. I’m wondering what the most robust approach to integrating complicated expressions might be. Consider the following:
```julia
import Symbolics
import SymbolicNumericIntegration
using SymbolicNumericIntegration: integrate
import QuadGK
Symbolics.@variables x
m = [x^n for n in range( 0, 4 )]
f = sin( π * x )^2 + cos( x ) - 1
integrand = m * f
```
```julia
julia> display( integrand )
5-element Vector{Symbolics.Num}:
-1 + cos(x) + sin(πx)^2
x*(-1 + cos(x) + sin(πx)^2)
(x^2)*(-1 + cos(x) + sin(πx)^2)
(x^3)*(-1 + cos(x) + sin(πx)^2)
(x^4)*(-1 + cos(x) + sin(πx)^2)
```
Using `integrate` from SymbolicsNumericIntegration.jl, which now supports definite integration, sometimes (often) returns the wrong answer:
```julia
julia> integrate( integrand[end], (x, 0, 1); detailed=false )
-0.06692331486014069
julia> integrate( integrand[end], (x, 0, 1); detailed=false )
-0.06692331486014069
julia> integrate( integrand[end], (x, 0, 1); detailed=false )
-0.009884419990295636
```
Cranking down the absolute tolerance doesn’t always help:
```julia
## abstol = 1e-12
julia> integrate( integrand[end], (x, 0, 1); detailed=false, abstol=1e-12 )
-0.06692331486014069
julia> integrate( integrand[end], (x, 0, 1); detailed=false, abstol=1e-12 )
-0.009884419990287524
julia> integrate( integrand[end], (x, 0, 1); detailed=false, abstol=1e-12 )
-0.009884419990296808
## abstol = 1e-20
julia> integrate( integrand[end], (x, 0, 1); detailed=false, abstol=1e-20 )
[ Info: Integration is partially successful. Pass `detailed = true` to `integrate` for details
```
This is especially annoying when broadcasting, as (seemingly) randomly some terms don’t integrate and thus return `nothing`:
```julia
var_tuple = fill( ( x, 0, 1 ), ( 5 ) )
julia> integrate.( integrand, var_tuple; detailed=false )
[ Info: Integration is partially successful. Pass `detailed = true` to `integrate` for details
[ Info: Integration is partially successful. Pass `detailed = true` to `integrate` for details
5-element Vector{Union{Nothing, Float64}}:
0.3414709848078965
0.13177329067603627
nothing
nothing
-0.06692331486014069
julia> integrate.( integrand, var_tuple; detailed=false )
[ Info: Integration is partially successful. Pass `detailed = true` to `integrate` for details
[ Info: Integration is partially successful. Pass `detailed = true` to `integrate` for details
5-element Vector{Union{Nothing, Float64}}:
0.3414709848078965
0.13177329067603627
nothing
-0.07826184164390249
nothing
```
I’ve generally found trying to tweak `kwargs` to be an exercise in futility and expression-dependent.
* * *
I’ve also tried using QuadGK.jl, resigning myself to using a purely numerical approach, but Symbolics don’t appear compatible with `QuadGK`?
```julia
julia> QuadGK.quadgk( integrand[end], 0, 1 )
ERROR: Sym (x^4)*(-1 + cos(x) + sin(πx)^2) is not callable. Use @syms (x^4)*(-1 + cos(x) + sin(πx)^2)(var1, var2,...) to create it as a callable.
```
It’s not clear to me how to create a `Symbolics` callable object in this case (where it’s desired to still be able to do things like `transpose` and `*` in other parts of the code):
```julia
Symbolics.@syms integrand_5( x )
julia> integrand_5(x) = x^2
ERROR: cannot define function integrand_5; it already has a value
julia> integrand_5
x^2
julia> integrand_5(2)
ERROR: Sym x^2 is not callable. Use @syms x^2(var1, var2,...) to create it as a callable.
```
Examples of `transpose` and `*` not working:
```julia
m_1(x) = x^0
m_2(x) = x^1
m_3(x) = x^2
m_4(x) = x^3
m_5(x) = x^4
M = [m_1, m_2, m_3, m_4, m_5]
julia> transpose( M )
1×5 transpose(::Vector{Function}) with eltype Union{}:
Error showing value of type LinearAlgebra.Transpose{Union{}, Vector{Function}}:
ERROR: MethodError: no method matching transpose(::typeof(m_1))
julia> M * f
ERROR: MethodError: no method matching *(::typeof(m_1), ::Symbolics.Num)
julia> reshape( M, 5, 1 ) * reshape( M, 1, 5 )
ERROR: MethodError: no method matching zero(::Type{Function})
```
* * *
Side note #1: the exact integral of `integrand[end]` can be found in SymPy / Matlab / Mathematica / etc., as:
-\frac{1}{10} + \frac{3}{4\pi^4} - \frac{1}{2\pi^2} + 13\sin(1) - 20\cos(1) \approx -0.009884419990295
* * *
Side note #2: the `Symbolics.jl` documentation claims that it supports “[Symbolic Linear Algebra](https://docs.sciml.ai/Symbolics/stable/#:~:text=Symbolic%20linear%20algebra%20(factorizations%2C%20inversion%2C%20determinants%2C%20eigencomputations%2C%20etc.))” however this rarely seems to work for me:
```julia
julia> var_tuple = fill( ( x, 0, 1 ), ( 5, 5 ) )
julia> A = integrate.( m * transpose( m ), var_tuple; detailed=false )
julia> var_tuple = fill( ( x, 0, 1 ), ( 5 ) )
julia> b = integrate.( m, var_tuple; detailed=false )
julia> A \ b
ERROR: MethodError: no method matching Float64(::Symbolics.Num)
```
I guess I expected Julia to be a bit more “composable” ala Matlab, isn’t that what multiple dispatch is supposed to provide?
Also, if you use matrix inversion the precision is a bit disappointing:
```julia
julia> inv( A ) * b
5-element Vector{Any}:
1.0000000000005969
-5.4569682106375694e-12
-1.8189894035458565e-12
3.637978807091713e-12
-3.637978807091713e-12
```
I thought the benefit of `Symbolic.Num` was exactness (to the extent possible)? That’s at least a major benefit of other symbolic engines. Consider Matlab:
```matlab
>> x = sym( "x", "real" )
>> m(x) = [x^0; x^1; x^2; x^3; x^4]
>> A = int( m * transpose( m ), [0, 1] )
>> b = int( m, [0, 1] )
>> sol = inv( A ) * b
sol =
1
0
0
0
0
>> isAlways( all( sol == sym( [1; 0; 0; 0; 0] ) ) )
ans =
logical
1
```
---
### Author:  [@TS-CUBED](https://discourse.julialang.org/u/TS-CUBED)
#### Post date: [October 17, 2025, 12:45pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/11 "2025-10-17T12:45:21Z")
If I shouldn’t revive old threads and you’d rather have me start a new one, please let me know. But I thought it might help others who stumble across this one.
Is SymbolicIntegration.jl the intended integration engine for Symbolics.jl? I played around with it and it seems to work nicely. However it limits Symbolics.jl to version 6.2, which doesn’t have some other features I need (symbolic\_solve for non-linear eqns).
Is it just a matter of giving it time and being patient with SymbolicIntegration.jl?
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [October 17, 2025, 12:48pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/12 "2025-10-17T12:48:54Z")
> [@TS-CUBED](#):
>
> Is [SymbolicIntegration.jl](https://juliaregistries.github.io/General/packages/redirect_to_repo/SymbolicIntegration) the intended integration engine for [Symbolics.jl](https://juliaregistries.github.io/General/packages/redirect_to_repo/Symbolics)?
Yes, this is now documented and released:
> **[Symbolic Integrals · Symbolics.jl](https://docs.sciml.ai/Symbolics/stable/manual/integration/)**
>
> Documentation for Symbolics.jl.
> **[SciML: Open Source Software for Scientific Machine Learning](https://sciml.ai/news/2025/10/10/SymbolicIntegration/)**
>
> Open Source Software for Scientific Machine Learning
> [@TS-CUBED](#):
>
> However it limits [Symbolics.jl](https://juliaregistries.github.io/General/packages/redirect_to_repo/Symbolics) to version 6.2, which doesn’t have some other features I need (symbolic\_solve for non-linear eqns).
It doesn’t have an upper bound, so in theory it shouldn’t be an issue. But I just checked and I think it has to do with the Nemo and AbstractAlgebra versions, so I’ll try and fix that in a little bit.
---
### Author:  [@jcbritobr](https://discourse.julialang.org/u/jcbritobr)
#### Post date: [November 5, 2025, 5:46pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/13 "2025-11-05T17:46:30Z")
> [@ChrisRackauckas](#):
>
> Implement a symbolic integration or convrt it to Sympy until it is implemented.
Hello, good afternoon.
We can now evaluate the results of Symbolics.Integral()?
Im trying to migrate from sympy, but cant find it in documentation.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [November 5, 2025, 11:33pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/14 "2025-11-05T23:33:58Z")
> **[Symbolic Integrals · Symbolics.jl](https://docs.sciml.ai/Symbolics/stable/manual/integration/)**
>
> Documentation for Symbolics.jl.
---
### Author:  [@jcbritobr](https://discourse.julialang.org/u/jcbritobr)
#### Post date: [November 9, 2025, 2:49pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/15 "2025-11-09T14:49:26Z")
Already read the docs. I can use SymbolicIntegration.jl, but Im refering to Symbolics.Integral, thats give me a symbol without evaluation. Anyway, im just using SymbolicIntegration.jl.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [November 9, 2025, 2:49pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/16 "2025-11-09T14:49:54Z")
The other is lazy
---
### Author:  [@jcbritobr](https://discourse.julialang.org/u/jcbritobr)
#### Post date: [November 9, 2025, 3:20pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/17 "2025-11-09T15:20:28Z")
Yes, but it cant be expanded like Differential. Does it?
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [November 9, 2025, 3:28pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/18 "2025-11-09T15:28:01Z")
It can’t, we should add that
---
### Author:  [@jcbritobr](https://discourse.julialang.org/u/jcbritobr)
#### Post date: [November 9, 2025, 3:29pm UTC](https://discourse.julialang.org/t/integrating-symbolic-functions/102545/19 "2025-11-09T15:29:33Z")
Thank you. That was giving-me some doubts.