quo_rem gives wrong answers for LaurentPolynomial_mpair

[Etn40ff]

The method quo_rem of LaurentPolynomial_mpair gives wrong answers:

sage: R.<x,y> = LaurentPolynomialRing(QQ)
sage: q, r = (1/x).quo_rem(y) ; q, r
(0, 1)
sage: q*y + r == 1/x
False

As remarked on sage-devel it is not clear which should be the preferred answer in this case: both (0, 1/x) and (1/(x*y), 0) are correct because y is a unit.

CC: @DaveWitteMorris

Component: algebra

Keywords: Laurent polynomials

Author: Dave Morris

Branch/Commit: 0b0db6d

Reviewer: Salvatore Stella

Issue created by migration from https://trac.sagemath.org/ticket/31257

[slel]

comment:1

sage: q, r = (1/x).quo_rem(y) ; q, r

Both (0, 1/x) and (1/(x*y), 0) are correct because y is a unit.

I would go for (0, 1/x).

[Etn40ff]

comment:2

Replying to @slel:

sage: q, r = (1/x).quo_rem(y) ; q, r

Both (0, 1/x) and (1/(x*y), 0) are correct because y is a unit.

I would go for (0, 1/x).

Integers and polynomials seem to adopt the other option:

sage: 10.quo_rem(-1)
(-10, 0)
sage: R.<x,y> = PolynomialRing(ZZ)
sage: (x+1).quo_rem(-1)
(-x - 1, 0)

[DaveWitteMorris]

Branch: public/31257

[DaveWitteMorris]

Author: Dave Morris

[DaveWitteMorris]

comment:4

The PR fixes two bugs.

A Laurent polynomial is stored as a pair: a polynomial, and a tuple of offsets for the exponents. The quo_rem method divides these two multivariable polynomials, but did not apply the appropriate offset to the remainder (although it did apply the correct offset to the quotient). This means that the remainder was usually wrong (not just in corner cases). Indeed, one of the doctests was incorrect:

sage: R.<s, t> = LaurentPolynomialRing(QQ)
sage: (s^-2-t^2).quo_rem(s-t)
(s + t, -s^4 + 1)

This is mathematically wrong, because the remainder should be s^-2 - s^2.

The other bug is less serious. It was caused by the non-uniqueness of the representation of a Laurent polynomial. For example, x^2 could be represented as x^2 with an offset of 0 in the exponent, or as the constant 1 with an offset of 2 in the exponent of x. Changing the representation will often change the result of quo_rem, so it could easily happen that g == g_1, but f.quo_rem(g) is not equal to f.quo_rem(g1). To fix this, I applied the normalize method in order to put self and right into standard form before applying the polynomial division. This should make the output of quo_rem unique.

A question, because I am not a cython expert (or python expert). The quo_rem method should not change self, so, instead of self.normalize(), I wrote

cdef LaurentPolynomial_mpair selfl = <LaurentPolynomial_mpair> self
selfl._normalize()

Is it true that this will normalize a copy of self, and leave self unchanged?

---

New commits:

0debcaa

trac 31257 fix quo_rem

[DaveWitteMorris]

Commit: 0debcaa

[Etn40ff]

comment:5

LGTM

In my understanding you are right: this should leave self unchanged. I am not sure this is the fastest solution though.

[Etn40ff]

Reviewer: Salvatore Stella

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1 and set ticket back to needs_review. New commits:

0b0db6d

make copies of input

[sagetrac-git]

Changed commit from 0debcaa to 0b0db6d

[Etn40ff]

comment:7

The code is now more in line with what happens in other places in the same file.

Talking of which, the univariate version of quo_rem does not normalize self before performing computations. Is this normalization required to get consistent answers?

[DaveWitteMorris]

comment:8

In the univariate case, Laurent polynomials are always stored in normalized form: "A Laurent series is a pair (u(t), n), where either u = 0 (to some precision) or u is a unit." The requirement that u is a unit means that the constant term is not 0, which is what it means to be normalized in the univariate case. Since the input is already normalized, there is no need to normalize it again.

But multivariate Laurent polynomials are not automatically normalized when they are created (or at other times).

[Etn40ff]

comment:9

Then I am fine with the fix.

[DaveWitteMorris]

comment:10

Thanks for the review!

For some reason, my explanation about the commit in comment:6 did not get uploaded, so I will paraphrase it now, for the record:

My original commit did change self. (Casting a pointer does not make a copy.) So it is necessary to explicitly make a copy of self and right.

[vbraun]

Changed branch from public/31257 to 0b0db6d