Symbolics.jl simplify function gives overflow error

I am trying out the Symbolics.jl package for the first time to avoid solving this 2x2 system by hand. However, I am getting an overflow error when trying to simplify the answer.

  1. Are there any strategies to avoid overflow? Use BigInt somehow or arrange as a single matrix equation maybe?
  2. Is there any way to prevent expressions from automatically simplifying upon input? I would like to check render(latexify(moment)) against my handwritten equation for input mistakes before it simplifies.
  3. Can I save a_solution and b_solution to a file for use later without resolving?
  4. Why isn’t solve_for exported from the Symbolics package?
using Symbolics
@variables a b P_o P_i r_o r_i
force = P_i*r_i^2 - P_o*r_o^2 + 2//3*a*(r_i^3-r_o^3) + b*(r_i^2-r_o^2)
moment = P_i*r_i^3 - P_o*r_o^3 + 3//4*a*(r_i^4-r_o^4) + b*(r_i^3-r_o^3)
(a_solution, b_solution) =
	Symbolics.solve_for(
	[0 ~ force, 0 ~ moment],
	[a, b],
	simplify = true # Solutions are found without this line.
	)
Error Message
OverflowError: -26873856 * -11145125560320 overflowed for type Int64
throw_overflowerr_binaryop(::Symbol, ::Int64, ::Int64)@checked.jl:154
checked_mul@checked.jl:288[inlined]
*(::Rational{Int64}, ::Rational{Int64})@rational.jl:334
_mul(::Type{Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}})@mult.jl:109
*@mult.jl:119[inlined]
_broadcast_getindex_evalf@broadcast.jl:670[inlined]
_broadcast_getindex@broadcast.jl:643[inlined]
getindex@broadcast.jl:597[inlined]
macro expansion@broadcast.jl:961[inlined]
macro expansion@simdloop.jl:77[inlined]
copyto!@broadcast.jl:960[inlined]
copyto!@broadcast.jl:913[inlined]
copy@broadcast.jl:885[inlined]
materialize(::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1}, Nothing, typeof(*), Tuple{Base.RefValue{DynamicPolynomials.Polynomial{true, Rational{Int64}}}, Vector{DynamicPolynomials.Polynomial{true, Rational{Int64}}}}})@broadcast.jl:860
_term_poly_mult(::DynamicPolynomials.Term{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::Function)@mult.jl:78
*@mult.jl:82[inlined]
_pseudo_rem(::MultivariatePolynomials.UniqueFactorizationDomain, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@division.jl:120
pseudo_rem(::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@division.jl:105
primitive_univariate_gcd(::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:390
univariate_gcd(::MultivariatePolynomials.UniqueFactorizationDomain, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::DynamicPolynomials.Polynomial{true, DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:470
univariate_gcd@gcd.jl:465[inlined]
multivariate_gcd@gcd.jl:332[inlined]
deflated_gcd(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:237
gcd(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:128
_simplifier(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:500
(::MultivariatePolynomials.var"#106#107"{MultivariatePolynomials.GeneralizedEuclideanAlgorithm})(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}})@gcd.jl:528
(::Base.BottomRF{MultivariatePolynomials.var"#106#107"{MultivariatePolynomials.GeneralizedEuclideanAlgorithm}})(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}})@reduce.jl:81
_foldl_impl(::Base.BottomRF{MultivariatePolynomials.var"#106#107"{MultivariatePolynomials.GeneralizedEuclideanAlgorithm}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::Vector{DynamicPolynomials.Polynomial{true, Rational{Int64}}})@reduce.jl:62
foldl_impl(::Base.BottomRF{MultivariatePolynomials.var"#106#107"{MultivariatePolynomials.GeneralizedEuclideanAlgorithm}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::Vector{DynamicPolynomials.Polynomial{true, Rational{Int64}}})@reduce.jl:48
mapfoldl_impl(::typeof(identity), ::MultivariatePolynomials.var"#106#107"{MultivariatePolynomials.GeneralizedEuclideanAlgorithm}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::Vector{DynamicPolynomials.Polynomial{true, Rational{Int64}}})@reduce.jl:44
_mapreduce_dim(::Function, ::Function, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::Vector{DynamicPolynomials.Polynomial{true, Rational{Int64}}}, ::Colon)@reducedim.jl:327
#mapreduce#725@reducedim.jl:322[inlined]
#reduce#727@reducedim.jl:371[inlined]
content@gcd.jl:527[inlined]
deflated_gcd(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:231
gcd(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::MultivariatePolynomials.GeneralizedEuclideanAlgorithm)@gcd.jl:128
gcd@gcd.jl:112[inlined]
_gcd(::DynamicPolynomials.Polynomial{true, Rational{Int64}}, ::DynamicPolynomials.Polynomial{true, Rational{Int64}})@polyform.jl:373
gcd(::SymbolicUtils.PolyForm{Real, Nothing}, ::SymbolicUtils.PolyForm{Real, Nothing})@polyform.jl:89
_gcd(::SymbolicUtils.PolyForm{Real, Nothing}, ::SymbolicUtils.PolyForm{Real, Nothing})@polyform.jl:373
rm_gcds(::Vector{SymbolicUtils.PolyForm{Real, Nothing}}, ::Vector{SymbolicUtils.PolyForm{Real, Nothing}})@polyform.jl:486
simplify_div(::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@polyform.jl:266
(::SymbolicUtils.var"#sdiv#115")(::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@polyform.jl:307
var"#_#83"(::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, ::ComposedFunction{SymbolicUtils.var"#sdiv#115", typeof(SymbolicUtils.quick_cancel)}, ::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@operators.jl:1085
(::ComposedFunction{SymbolicUtils.var"#sdiv#115", typeof(SymbolicUtils.quick_cancel)})(::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@operators.jl:1085
(::Metatheory.Rewriters.Walk{:post, ComposedFunction{SymbolicUtils.var"#sdiv#115", typeof(SymbolicUtils.quick_cancel)}, typeof(SymbolicUtils.frac_similarterm), false})(::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@Rewriters.jl:200
var"#simplify_fractions#114"(::Bool, ::typeof(SymbolicUtils.simplify_fractions), ::SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing})@polyform.jl:309
simplify_fractions@polyform.jl:303[inlined]
_broadcast_getindex_evalf@broadcast.jl:670[inlined]
_broadcast_getindex@broadcast.jl:643[inlined]
_getindex@broadcast.jl:667[inlined]
_broadcast_getindex@broadcast.jl:642[inlined]
getindex@broadcast.jl:597[inlined]
copyto_nonleaf!(::Vector{SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing}}, ::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1}, Tuple{Base.OneTo{Int64}}, typeof(SymbolicUtils.simplify), Tuple{Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1}, Nothing, typeof(SymbolicUtils.simplify_fractions), Tuple{Base.Broadcast.Extruded{Vector{SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing}}, Tuple{Bool}, Tuple{Int64}}}}}}, ::Base.OneTo{Int64}, ::Int64, ::Int64)@broadcast.jl:1055
copy@broadcast.jl:907[inlined]
materialize(::Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1}, Nothing, typeof(SymbolicUtils.simplify), Tuple{Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1}, Nothing, typeof(SymbolicUtils.simplify_fractions), Tuple{Vector{SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing}}}}}})@broadcast.jl:860
var"#solve_for#242"(::Bool, ::Bool, ::typeof(Symbolics.solve_for), ::Vector{Symbolics.Equation}, ::Vector{Symbolics.Num})@linear_algebra.jl:110
top-level scope@Local: 1[inlined]

Here is the source of the equations I am solving in case that matters.


Do you need the rational // form for the solution? If you evaluate it, you loose the infinite precision.
That said, the code runs with the simplify line if I change // to normal division, or mulitiply by the divisors.
The error message hints to come from problems (cancelations?) involving gcd. I would suggest raising an issue.

1 Like

Without looking into the details I tried

using Symbolics
@variables a b P_o P_i r_o r_i
force = P_i*r_i^2 - P_o*r_o^2 + BigInt(2)//3*a*(r_i^3-r_o^3) + b*(r_i^2-r_o^2)
moment = P_i*r_i^3 - P_o*r_o^3 + BigInt(3)//4*a*(r_i^4-r_o^4) + b*(r_i^3-r_o^3)
(a_solution, b_solution) =
	Symbolics.solve_for(
	[0 ~ force, 0 ~ moment],
	[a, b],
	simplify = true # Solutions are found without this line.
	)

and got

2-element Vector{SymbolicUtils.Div{Real, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, SymbolicUtils.Add{Real, Int64, Dict{Any, Number}, Nothing}, Nothing}}:
 (P_i*(r_i^2)*(r_o^2) - P_o*(r_i^2)*(r_o^2)) / ((1//12)*(r_i^5) + (1//12)*r_o*(r_i^4) + (2//3)*(r_i^2)*(r_o^3) - (1//12)*(r_o^5) - (2//3)*(r_i^3)*(r_o^2) - (1//12)*r_i*(r_o^4))
 ((72//1)*P_i*(r_i^5) + (72//1)*P_i*r_o*(r_i^4) + (648//1)*P_i*(r_i^2)*(r_o^3) + (72//1)*P_i*(r_i^3)*(r_o^2) - (72//1)*P_o*(r_o^5) - (72//1)*P_o*r_i*(r_o^4) - (648//1)*P_o*(r_i^3)*(r_o^2) - (72//1)*P_o*(r_i^2)*(r_o^3)) / ((72//1)*(r_o^5) + (576//1)*(r_i^3)*(r_o^2) + (72//1)*r_i*(r_o^4) - (72//1)*(r_i^5) - (72//1)*r_o*(r_i^4) - (576//1)*(r_i^2)*(r_o^3))
1 Like