Consider changing backend for `PolyForm` arithmetic type promotion

[bowenszhu]

Short description

The current implementation of PolyForm division relies on DynamicPolynomials.Polynomial, which in turn utilizes MutableArithmetics.jl for type promotion. MutableArithmetics.jl, designed for numerical software such as JuMP, promotes division of two integers to a Float64. This behavior might not be ideal for symbolic computations where preserving integer types is often desired.

Detailed description

If the input expression contains a / operation, simplify and simplify_fractions call simplify_div. https://github.com/JuliaSymbolics/SymbolicUtils.jl/blob/6713fa0cf539c611b97aa6bed9e0fc58923a4f8d/src/simplify.jl#L16-L43 https://github.com/JuliaSymbolics/SymbolicUtils.jl/blob/6713fa0cf539c611b97aa6bed9e0fc58923a4f8d/src/polyform.jl#L327-L337

simplify_div transforms both numerator and denominator to PolyForms, and removes their greatest common divisor. https://github.com/JuliaSymbolics/SymbolicUtils.jl/blob/6713fa0cf539c611b97aa6bed9e0fc58923a4f8d/src/polyform.jl#L277-L280

rm_gcds involves the division of PolyForms. https://github.com/JuliaSymbolics/SymbolicUtils.jl/blob/6713fa0cf539c611b97aa6bed9e0fc58923a4f8d/src/polyform.jl#L534-L543 https://github.com/JuliaSymbolics/SymbolicUtils.jl/blob/6713fa0cf539c611b97aa6bed9e0fc58923a4f8d/src/polyform.jl#L86 As shown in the above code, our current implementation of PolyForm division relies on DynamicPolynomials.Polynomial division via div(x.p, y.p).

JuliaAlgebra/MultivariatePolynomials.jl uses jump-dev/MutableArithmetics.jl for arithmetic type promotion. https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/blob/9eb6bf4973e27fdbf6f53d0d2043b476ad2a38c5/src/division.jl#L148-L150 https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/blob/9eb6bf4973e27fdbf6f53d0d2043b476ad2a38c5/src/division.jl#L388-L390 https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/blob/9eb6bf4973e27fdbf6f53d0d2043b476ad2a38c5/src/division.jl#L379-L384 https://github.com/jump-dev/MutableArithmetics.jl/blob/7344c9468072bd663feec7e4cdc1a0e98e34ea19/src/interface.jl#L112-L114 https://github.com/jump-dev/MutableArithmetics.jl/blob/7344c9468072bd663feec7e4cdc1a0e98e34ea19/src/interface.jl#L38-L44

In the last piece of code above, for inferring the resulting type of two Ints with the / operation, it obtains Float64 instead of Int or Rational by essentially executing typeof(Int(0) / Int(1)). This is reasonable for numerical software, but we might want it to behave differently in our symbolic-focused software.

[hersle]

The simple change in JuliaAlgebra/MultivariatePolynomials.jl#294 would fix all examples in JuliaSymbolics/Symbolics.jl#1083. My knowledge of that library is limited, but I think it might break other things there, and not be ready to merge (?)

[blegat]

I didn't want to force div of polynomials with integer coefficients to be floating point, it was just a case I didn't give much thought on. Indeed, div and rem for polynomials only make sense when the coefficient type is a field. For instance divrem(5x, 3x) should be 5/3 with zero remainder, not 1 with remainder 2x because the leading monomial of the remainder cannot be divisible by the leading monomial of the divisor. There is also pseudo_rem in case you have non-field coefficient like Int. For that reason, I assumed the user that needs to divide polynomials with integer coefficients would promote it to either rational or float in case there was a strong opinion on what the result should be. So there is a quick fix on your end:

julia> function to_field(p)
           T = coefficient_type(p)
           if T <: Integer
               return map_coefficients(c -> Rational{T}(c), p)
           else
               return p
           end
       end
to_field (generic function with 1 method)

julia> to_field(2x + 1)
1//1 + 2//1x

Now, defaulting to Rational could make sense, I'm having an attempt in https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/pull/296