avoiding integer overflow in symbolic calculations?

[stevengj]

This is an example from our matrix calculus class that worked in SymPy but is failing with Symbolics:

function Babylonian(x; N = 10) 
    t = (1+x)/2
    for i = 2:N; t=(t + x/t)/2  end    
    t
end  

using Symbolics
@variables x
simplify(Babylonian(x,N=4))

which gives

julia> simplify(Babylonian(x,N=4))
ERROR: OverflowError: -11118795390271 * 4572599641 overflowed for type Int64

followed by a long stacktrace.

Is there a way to tell Symbolics that I want it to use BigInt arithmetic here? e.g. shouldn't (1+x)/2 use BigInt if x is symbolic?

[ChrisRackauckas]

It just follows Julia's rules, so (1+x)/2 uses regular integers but (1+x)/big(2) uses BigInt.

function Babylonian(x; N = 10) 
    t = (1+x)/2
    for i = 2:N; t=(t + x/t)/big(2)  end    
    t
end  

using Symbolics
@variables x
Symbolics.derivative(simplify(Babylonian(x,N=4)),x)
simplify(Symbolics.derivative(simplify(Babylonian(x,N=4)),x))

julia> Symbolics.derivative(simplify(Babylonian(x,N=4)),x)
(15 + 454x + 3003(x^2) + 6432(x^3) + 5005(x^4) + 1362(x^5) + 105(x^6)) / (128(1 + x)*((1//4) + (3//2)*x + (1//4)*(x^2))*((1//16) + (7//4)*x + (35//8)*(x^2) + (7//4)*(x^3) + (1//16)*(x^4))) - (128((1//4) + (3//2)*x + (1//4)*(x^2))*((1//16) + (7//4)*x + (35//8)*(x^2) + (7//4)*(x^3) + (1//16)*(x^4)) + 128(1 + x)*((1//4) + (3//2)*x + (1//4)*(x^2))*((7//4) + (35//4)*x + (21//4)*(x^2) + (1//4)*(x^3)) + 128(1 + x)*((3//2) + (1//2)*x)*((1//16) + (7//4)*x + (35//8)*(x^2) + (7//4)*(x^3) + (1//16)*(x^4)))*((15x + 227(x^2) + 1001(x^3) + 1608(x^4) + 1001(x^5) + 227(x^6) + 15(x^7)) / (16384((1 + x)^2)*(((1//4) + (3//2)*x + (1//4)*(x^2))^2)*(((1//16) + (7//4)*x + (35//8)*(x^2) + (7//4)*(x^3) + (1//16)*(x^4))^2)))

julia> simplify(Symbolics.derivative(simplify(Babylonian(x,N=4)),x))
(960.0 + 29056.0x + 438592.0(x^2) + 3.523328e+06(x^3) + 1.6168832e+07(x^4) + 4.342272e+07(x^5) + 7.0727424e+07(x^6) + 7.0659072e+07(x^7) + 4.3384384e+07(x^8) + 1.6060544e+07(x^9) + 3.474368e+06(x^10) + 396032.0(x^11) + 19072.0(x^12)) / (524288((1 + x)^2)*(((1//4) + (3//2)*x + (1//4)*(x^2))^2)*(((1//16) + (7//4)*x + (35//8)*(x^2) + (7//4)*(x^3) + (1//16)*(x^4))^2))

works, though the fact that the last simplify changes to floats is a bug. @shashi do you know why that might be happening?

[stevengj]

Julia evaluation for (1+x)/2 behaves differently depending on the type of x, so I'm not sure what you mean by "It just follows Julia's rules". I don't want to use (1+x)/BigInt(2) because if x is a Float64 I want to keep it in that precision instead of promoting to BigFloat. You are forcing me to write a different Babylonian method just for Symbolics.

What I would like would be a promotion rule (maybe enabled by a flag or a type of symbolic variable?), that allows one to use (1+x)/2, and says that when you combine integers with a symbolic expression you want to keep everything exact (i.e. promote to BigInt, use rationals, etcetera), ala Mathematica or SymPy.

Is there a way to attach a type to a symbolic variable? e.g. @variables x::Rational{BigInt} would make it do promotion as if x were an unknown value of that type, keeping all coefficients exact. The current behavior seems more like @variables x::Float64 (float coefficients) or @variables x::Rational{Int} (overflow).