# Test of expand and substitute in Symbolics.jl, SymEngine.jl and Mathematica
**URL:** https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470
**Category:** Performance
**Created:** [March 6, 2022, 10:18am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470 "2022-03-06T10:18:21Z")
**Posts on this page:** 19
**Page:** 1
### Author:  [@gangchern](https://discourse.julialang.org/u/gangchern)
#### Post date: [March 6, 2022, 10:18am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/1 "2022-03-06T10:18:21Z")
My computer is Mac Pro M1, 2020.
General Comments: for expand, SymEngine is nearly as good as Mathematica, Symbolics‘s “expand” function still needs to be greatly improved. substitute are similar for the three CAS in performance.
1. expand,
Symbolics (code+time)
```julia
f=(x+big(2)y)^6000|>expand; # 64.2s
f=(x+big(2)y)^100000|>expand; # looks infinity
```
SymEngine
```julia
f=(x+2y)^6000|>expand; #0.5s
f=(x+2y)^100000|>expand; # 10.6s
```
Mathematica takes 0.01s and 3.1s respectively.
1. substitute
Symbolics takes 0.9s
```julia
substitute(f,Dict(x=>big(1),y=>big(2)));
```
SymEngine takes 0.4s
```julia
subs(f,Dict(x=>1,y=>2));
```
Mathematica takes 0.06s
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [March 6, 2022, 12:42pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/2 "2022-03-06T12:42:30Z")
Yes, it does need to be improved. [https://github.com/JuliaSymbolics/SymbolicUtils.jl/pull/432](https://github.com/JuliaSymbolics/SymbolicUtils.jl/pull/432) should help a lot but is not released right now to Symbolics.
---
### Author:  [@gangchern](https://discourse.julialang.org/u/gangchern)
#### Post date: [March 9, 2022, 5:29am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/3 "2022-03-09T05:29:24Z")
great news! Looking forward to that.
---
### Author:  [@gangchern](https://discourse.julialang.org/u/gangchern)
#### Post date: [June 5, 2022, 5:41am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/4 "2022-06-05T05:41:16Z")
Now I tested again. Symbolics is much faster than SymEngine now
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 5, 2022, 7:21pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/5 "2022-06-05T19:21:12Z")
```julia
julia> using Symbolics
julia> @variables x, y
julia> @time f = (x + y)^6000 |>expand;
14.399519 seconds (30.27 M allocations: 3.254 GiB, 3.97% gc time, 25.09% compilation time)
julia> using SymEngine
julia> @vars x, y
julia> f = (x + y)^6000 |> expand;
julia> @time f = (x + y)^6000 |>expand;
0.012988 seconds (36.10 k allocations: 12.602 MiB)
```
In my Mac M1, `SymEngine` is faster than `Symbolics`.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [June 5, 2022, 10:56pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/6 "2022-06-05T22:56:17Z")
Is that using the unityper form on master?
---
### Author:  [@gangchern](https://discourse.julialang.org/u/gangchern)
#### Post date: [June 6, 2022, 3:49am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/7 "2022-06-06T03:49:43Z")
```julia
using Symbolics
@variables x, y
@time f = (x + 2y)^50000 |>expand;
0.026855 seconds (70.92 k allocations: 6.684 MiB, 30.24% gc time, 9.12% compilation time)
```
Strange! looks too good to be true! Maybe I did a wrong test.
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:31pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/8 "2022-06-08T20:31:32Z")
```julia
julia> @time f = (x + 2y)^50000 |>expand;
0.022756 seconds (51.71 k allocations: 4.784 MiB)
julia> @time f = (x + y)^6000 |>expand;
10.740392 seconds (15.51 M allocations: 2.443 GiB, 4.71% gc time)
```
It’s so strange.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [June 8, 2022, 8:36pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/9 "2022-06-08T20:36:29Z")
That’s because of the type-instability that unityper fixes.
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:38pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/10 "2022-06-08T20:38:47Z")
why is `(x + 2y)^50000 |>expand` is faster than `f = (x + y)^50000`?
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:41pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/11 "2022-06-08T20:41:27Z")
```julia
julia> @time f = (x + 2y)^50000 |>expand;
0.023723 seconds (51.76 k allocations: 4.787 MiB)
julia> @time f = (x + y)^50000 |>expand;
```
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:51pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/12 "2022-06-08T20:51:41Z")
```julia
julia> @time f = (x + 4y)^50000 |>expand;
0.005870 seconds (16.61 k allocations: 1.465 MiB)
julia> @time f = (x + 6y)^50000 |>expand;
0.044559 seconds (109.41 k allocations: 8.064 MiB, 56.32% compilation time)
julia> @time f = (x + 8y)^50000 |>expand;
0.051483 seconds (40.72 k allocations: 2.619 MiB, 92.70% compilation time)
julia> @time f = (x + 10y)^50000 |>expand;
0.028289 seconds (51.71 k allocations: 4.784 MiB, 31.20% gc time)
```
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [June 8, 2022, 8:52pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/13 "2022-06-08T20:52:41Z")
You may have hit a compilation phase?
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:54pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/14 "2022-06-08T20:54:37Z")
It seems that when the coefficient of `y` is even, it will be very fast.
But when the coefficient is odd, it will run for a long long long time.
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 8:57pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/15 "2022-06-08T20:57:00Z")
```julia
julia> for i in 1:10
@time f = (x + i * y)^5000 |>expand
end
5.883341 seconds (10.36 M allocations: 1.695 GiB, 6.25% gc time, 0.07% compilation time)
0.008549 seconds (39.72 k allocations: 3.583 MiB)
5.866237 seconds (10.36 M allocations: 1.694 GiB, 5.82% gc time)
0.002496 seconds (13.64 k allocations: 1.154 MiB)
6.197537 seconds (10.36 M allocations: 1.694 GiB, 9.01% gc time)
0.008952 seconds (39.72 k allocations: 3.582 MiB)
5.930975 seconds (10.36 M allocations: 1.694 GiB, 5.82% gc time)
0.001464 seconds (7.69 k allocations: 619.297 KiB)
6.237053 seconds (10.36 M allocations: 1.694 GiB, 10.19% gc time)
0.011223 seconds (39.72 k allocations: 3.582 MiB)
```
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 8, 2022, 9:05pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/16 "2022-06-08T21:05:02Z")
```julia
julia> for i in 1:10
@time f = (x + i * y)^10000 |>expand
end
38.717680 seconds (39.79 M allocations: 5.797 GiB, 26.24% gc time, 0.01% compilation time)
0.009729 seconds (42.83 k allocations: 3.898 MiB)
31.223451 seconds (39.79 M allocations: 5.797 GiB, 11.21% gc time)
0.002769 seconds (14.51 k allocations: 1.247 MiB)
32.939406 seconds (39.79 M allocations: 5.797 GiB, 13.77% gc time)
0.009396 seconds (42.83 k allocations: 3.898 MiB)
31.781937 seconds (39.79 M allocations: 5.797 GiB, 10.80% gc time)
0.043063 seconds (37.99 k allocations: 2.418 MiB, 95.46% compilation time)
32.121851 seconds (39.79 M allocations: 5.797 GiB, 10.63% gc time)
0.009469 seconds (42.83 k allocations: 3.898 MiB)
```
---
### Author:  [@contradict](https://discourse.julialang.org/u/contradict)
#### Post date: [June 8, 2022, 11:26pm UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/17 "2022-06-08T23:26:11Z")
The pattern goes away if you use floating point numbers. This looks like integer overflow to me.
```julia
julia> for i in 1.0:10.0
@time f=(x+i*y)^5000 |> expand
end
8.599273 seconds (10.36 M allocations: 1.696 GiB, 8.55% gc time)
8.491495 seconds (10.36 M allocations: 1.696 GiB, 9.21% gc time)
8.400282 seconds (10.36 M allocations: 1.696 GiB, 8.55% gc time)
8.377576 seconds (10.36 M allocations: 1.696 GiB, 8.16% gc time)
8.375058 seconds (10.37 M allocations: 1.696 GiB, 8.39% gc time, 0.09% compilation time)
8.483856 seconds (10.36 M allocations: 1.696 GiB, 8.57% gc time)
8.438680 seconds (10.36 M allocations: 1.696 GiB, 8.32% gc time)
8.378851 seconds (10.36 M allocations: 1.696 GiB, 8.54% gc time)
8.499839 seconds (10.36 M allocations: 1.696 GiB, 9.27% gc time)
8.608331 seconds (10.37 M allocations: 1.696 GiB, 8.55% gc time)
```
Interestingly,
```julia
julia> (2*y)^5000
ERROR: OverflowError: 4294967296 * 4294967296 overflowed for type Int64
```
but
```julia
julia> (x+2*y)^5000
(x + 2y)^5000
```
```julia
```
---
### Author:  [@Raymond](https://discourse.julialang.org/u/Raymond)
#### Post date: [June 9, 2022, 12:13am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/18 "2022-06-09T00:13:28Z")
```julia
julia> (x + y)^66 |> expand
x^66 + y^66 + 66y*(x^65) + 7219428434016265740(x^33)*(y^33) + 4472995859186094240(x^37)*(y^29) + 6406484391866534976(x^31)*(y^35) + 3413602103063071920(x^28)*(y^38) + 4922879481520(x^12)*(y^54) + 2515943049305400(x^49)*(y^17) + 77413632286320(x^52)*(y^14) + 17302625882942400(x^19)*(y^47) + 3413602103063071920(x^38)*(y^28) + 17302625882942400(x^47)*(y^19) + 37014131440(x^9)*(y^57) + 6406484391866534976(x^35)*(y^31) + 4922879481520(x^54)*(y^12) + 2450791253481179840(x^39)*(y^27) + 348524321356411200(x^23)*(y^43) + 77413632286320(x^14)*(y^52) + 40661170824914640(x^20)*(y^46) + 268367258592576(x^15)*(y^51) + 5743572120(x^8)*(y^58) + 20448884000160(x^53)*(y^13) + 45760(x^63)*(y^3) + 855420636763836(x^16)*(y^50) + 855420636763836(x^50)*(y^16) + 2145(x^64)*(y^2) + 210980549208(x^56)*(y^10) + 1049058207282797712(x^41)*(y^25) + 2145(x^2)*(y^64) + 1074082795968(x^55)*(y^11) + 5516694892996182896(x^30)*(y^36) + 90858768(x^6)*(y^60) + 1074082795968(x^11)*(y^55) + 210980549208(x^10)*(y^56) + 90858768(x^60)*(y^6) + 1049058207282797712(x^25)*(y^41) + 7007092303604022630(x^34)*(y^32) + 720720(x^62)*(y^4) + 6848956078664700(x^18)*(y^48) + 268367258592576(x^51)*(y^15) + 348524321356411200(x^43)*(y^23) + 778789440(x^59)*(y^7) + 7007092303604022630(x^32)*(y^34) + 4472995859186094240(x^29)*(y^37) + 37014131440(x^57)*(y^9) + 2515943049305400(x^17)*(y^49) + 45760(x^3)*(y^63) + 1654284096099796392(x^26)*(y^40) + 89067326568860640(x^21)*(y^45) + 182183167981760400(x^22)*(y^44) + 5743572120(x^58)*(y^8) + 89067326568860640(x^45)*(y^21) + 778789440(x^7)*(y^59) + 5516694892996182896(x^36)*(y^30) + 8936928(x^5)*(y^61) + 20448884000160(x^13)*(y^53) + 6848956078664700(x^48)*(y^18) + 720720(x^4)*(y^62) + 182183167981760400(x^44)*(y^22) + 2450791253481179840(x^27)*(y^39) + 66x*(y^65) + 40661170824914640(x^46)*(y^20) + 624439409096903400(x^24)*(y^42) + 624439409096903400(x^42)*(y^24) + 8936928(x^61)*(y^5) + 1654284096099796392(x^40)*(y^26)
julia> (x + y)^67 |> expand
x^67 + y^67 + 67y*(x^66) + 766480(x^4)*(y^63) + 766480(x^63)*(y^4) + 7886597962249166160(x^29)*(y^38) + 972963730453314600(x^43)*(y^24) + 3371363686069236(x^50)*(y^17) + 3371363686069236(x^17)*(y^50) + 129728497393775280(x^46)*(y^21) + 99795696(x^6)*(y^61) + 4105075349580976232(x^40)*(y^27) + 972963730453314600(x^24)*(y^43) + 67x*(y^66) + 2703342303382594104(x^41)*(y^26) + 1673497616379701112(x^42)*(y^25) + 5996962277488(x^12)*(y^55) + 99795696(x^61)*(y^6) + 2703342303382594104(x^26)*(y^41) + 1673497616379701112(x^25)*(y^42) + 869648208(x^7)*(y^60) + 2211(x^65)*(y^2) + 9364899127970100(x^49)*(y^18) + 345780890878896(x^52)*(y^15) + 247994680648(x^57)*(y^10) + 47905(x^3)*(y^64) + 9657648(x^62)*(y^5) + 271250494550621040(x^22)*(y^45) + 42757703560(x^9)*(y^58) + 2211(x^2)*(y^65) + 247994680648(x^10)*(y^57) + 42757703560(x^58)*(y^9) + 4105075349580976232(x^27)*(y^40) + 1285063345176(x^56)*(y^11) + 47905(x^64)*(y^3) + 97862516286480(x^14)*(y^53) + 57963796707857040(x^47)*(y^20) + 97862516286480(x^53)*(y^14) + 1123787895356412(x^16)*(y^51) + 24151581961607100(x^48)*(y^19) + 57963796707857040(x^20)*(y^47) + 1285063345176(x^11)*(y^56) + 5996962277488(x^55)*(y^12) + 530707489338171600(x^23)*(y^44) + 5864393356544251760(x^39)*(y^28) + 271250494550621040(x^45)*(y^22) + 869648208(x^60)*(y^7) + 6522361560(x^8)*(y^59) + 5864393356544251760(x^28)*(y^39) + 9364899127970100(x^18)*(y^49) + 25371763481680(x^54)*(y^13) + 25371763481680(x^13)*(y^54) + 129728497393775280(x^21)*(y^46) + 9657648(x^5)*(y^62) + 24151581961607100(x^19)*(y^48) + 345780890878896(x^15)*(y^52) + 530707489338171600(x^44)*(y^23) + 1123787895356412(x^51)*(y^16) + 7886597962249166160(x^38)*(y^29) + 6522361560(x^59)*(y^8) - 4220223336089263246(x^33)*(y^34) - 5033167378238994010(x^32)*(y^35) - 6523564788846833744(x^31)*(y^36) - 6523564788846833744(x^36)*(y^31) - 8457053321527274480(x^37)*(y^30) - 4220223336089263246(x^34)*(y^33) - 8457053321527274480(x^30)*(y^37) - 5033167378238994010(x^35)*(y^32)
```
it seems that when the exponent is larger than 66, the result is wrong.
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [June 9, 2022, 2:38am UTC](https://discourse.julialang.org/t/test-of-expand-and-substitute-in-symbolics-jl-symengine-jl-and-mathematica/77470/19 "2022-06-09T02:38:24Z")
Yeah it’s overflow. Change the integer type