Modulo with polynomials is inconsistent

[SoongNoonien]

Hi!
While playing around with some modulo computations I stumbled arccos some weird thinks in Nemo and AbstractAlgebra. Apparently Nemo interprets mod differently than AbstractAlgebra:

julia> using Nemo

Welcome to Nemo version 0.34.7

Nemo comes with absolutely no warranty whatsoever


julia> R,(x,y)=ZZ[:x,:y]
(Multivariate Polynomial Ring in x, y over Integer Ring, ZZMPolyRingElem[x, y])

julia> mod(R(-1),R(2))
-1

julia> rem(R(-1),R(2),RoundDown)
-1

While

julia> using AbstractAlgebra

julia> R,(x,y)=ZZ[:x,:y]
(Multivariate polynomial ring in 2 variables over integers, AbstractAlgebra.Generic.MPoly{BigInt}[x, y])

julia> mod(R(-1),R(2))
1

julia> rem(R(-1),R(2),RoundDown)
1

works as I would expect.

[fingolfin]

This happens in FLINT, specifically fmpz_mpoly_divrem, which says this:

.. function:: void fmpz_mpoly_divrem(fmpz_mpoly_t Q, fmpz_mpoly_t R, const fmpz_mpoly_t A, const fmpz_mpoly_t B, const fmpz_mpoly_ctx_t ctx)

Set Q and R to the quotient and remainder of A divided by B. The monomials in R divisible by the leading monomial of B will have coefficients reduced modulo the absolute value of the leading coefficient of B.

So it claims that "coefficients [are] reduced modulo [...] the leading coefficient of B."

Combine this with the observation that

julia> rem(ZZ(-1),ZZ(2))
-1

julia> rem(-1,2)
-1

So perhaps it reduces with RoundToZero, not RoundDown?

This is not quite trivial to change for us. Consider:

julia> R,(x,y) = ZZ[:x,:y];

julia> mod(-x^2-3x, 2x)
-x^2 - x

julia> mod(x^2+3x, 2x)
x^2 + x

julia> mod(x^2-3x, 2x)
x^2 - x

julia> mod(x^2+x, 2x)
x^2 + x

julia> mod(x^2-x, 2x)
x^2 - x

I'd expect each to give the same result, x^2+x.

It is even weirder (to me at least) in the univariate case, because there -3x mod 2x is reduced to x but -x mod 2x stays -x... Huh?!:

julia> R,x = ZZ[:x]
(Univariate polynomial ring in x over ZZ, x)

julia> mod(-x^2-3x, 2x)
-x^2 + x

julia> mod(x^2+3x, 2x)
x^2 + x

julia> mod(x^2-3x, 2x)
x^2 + x

julia> mod(x^2+x, 2x)
x^2 + x

julia> mod(x^2-x, 2x)
x^2 - x