If you want it exact, don’t ask for the Float64 result?
julia> ex = Term(sqrt,[big(2)])
sqrt(2)
julia> eval(toexpr(ex))
1.414213562373095048801688724209698078569671875376948073176679737990732478462102
It just follows Julia’s semantics.
I see. Yeah there’s a more interesting question to be answered here, which is what kind of number is the default symbol? The default is @variables x::Real
, specifically not ::AbstractFloat
. Since these rules are all type-based, I think simplifying x-x -> 0
on ::Real
is the semantics you’d expect, while x-x -> 0 * x
should be what happens on ::AbstractFloat
. This keeps two algebras distinct: real numbers obey commutative and associative laws, they are a field, etc. while floating point numbers have a different set of rules.
The reason for keeping this difference is that the set of transformations that keep exact floating point semantics is already in LLVM, so I don’t think there’s too much to gain there. But it would be a good opt-in: “make simplifications only do what’s allowed on floating point numbers” instead of “make simplifications only do what’s allowed on real numbers”. The latter seems to be what people normally want on models, or at least out of a CAS, though there are times when exact floating point semantics need to be followed (in which case we have @register
for someone to register a function and keep some calculation intact). What I want to do is add a call to Herbie, so that in many cases we aren’t just applying rules that change floating point behavior but actually improve it, so even further down the path of “keep floating point separate”.