sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.38k stars 470 forks source link

clean univariate Laurent polynomials #25524

Open videlec opened 6 years ago

videlec commented 6 years ago
  1. faster creation of new Laurent polynomial via a cdef _new_c method
  2. Make rich comparison compatible with polynomials
  3. make Laurent polynomial complient with coercion for / so that the result is always a rational fration.
  4. allow construction of Laurent polynomial from rational fraction

As a consequence of 2 and 3, many doctests had to be corrected in sage/algebras/.

CC: @tscrim @vinklein

Component: basic arithmetic

Author: Vincent Delecroix

Branch/Commit: u/vdelecroix/25524 @ 228529a

Reviewer: Travis Scrimshaw

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

videlec commented 6 years ago

New commits:

028a39e25524: cleaning Laurent polynomials
20291d825524: some fixes to sage/algebras
9cacb5a25524: fix doctests in algebras
videlec commented 6 years ago

Commit: 9cacb5a

videlec commented 6 years ago

Branch: u/vdelecroix/25524

tscrim commented 6 years ago
comment:2

I disagree with 1/x not being a Laurent polynomial (and I think by extension ~x). This is far too common for people to do and then have to deal with rational functions for something that looks like a Laurent polynomial. IIRC, there were some sage-devel discussions on this. If you do want to make this change, then you need to ask on sage-devel.

I would also appreciate it if you could put back in the linebreaks on the long outputs.

The _new_c also does not mitigate the subclass issue. In fact, it makes it worse if I want to do a Python subclass as I now have no way to override _new_c used in those other methods. I would instead use

cdef type t = type(self)
cdef LaurentPolynomial_univariate p = t.__new__(t)

following sage.groups.perm_gps.permgroup_element.PermutationGroupElement._new_c.

videlec commented 6 years ago
comment:3

Replying to @tscrim:

I disagree with 1/x not being a Laurent polynomial (and I think by extension ~x). This is far too common for people to do and then have to deal with rational functions for something that looks like a Laurent polynomial.

People are familiar with integers and

sage: parent(1 / 1)
Rational Field

Here 1 / 1 looks like an integer. Also people are more familiar with polynomials rather than Laurent polynomials and

sage: R = QQ['x']
sage: parent(R.one() / R.one())
Fraction Field of Univariate Polynomial Ring in x over Rational Field

Still R.one() / R.one() looks like an integer.

IIRC, there were some sage-devel discussions on this. If you do want to make this change, then you need to ask on sage-devel.

I am willing to ask on sage-devel as the current situation is a bug with respect to

The parent of `a / b` should only depend on the parent of `a` and `b`.

I can imagine two alternatives

  1. Allow a / b only when b is invertible which is what you seem to suggest. As a consequence 1 / (1 + q) will raise an error and one would have to explicitely convert to the fraction field to make it work.
  2. Make a / b always return a rational fraction. This is the behavior of all parents in Sage (or at least integers and polynomials). This is also what I proposed in the current branch.

I would also appreciate it if you could put back in the linebreaks on the long outputs.

This I can do. The changes in output look ok? They seem to come from the change in ordering for Laurent polynomials.

The _new_c also does not mitigate the subclass issue. In fact, it makes it worse if I want to do a Python subclass as I now have no way to override _new_c used in those other methods. I would instead use

cdef type t = type(self)
cdef LaurentPolynomial_univariate p = t.__new__(t)

following sage.groups.perm_gps.permgroup_element.PermutationGroupElement._new_c.

I have to investigate how Cython is optimizing this.

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 6 years ago

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

2c4184725524: cleaning Laurent polynomials
19ca66625524: some fixes to sage/algebras
a9b43eb25524: fix doctests in algebras
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 6 years ago

Changed commit from 9cacb5a to a9b43eb

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 6 years ago

Branch pushed to git repo; I updated commit sha1. New commits:

228529a25524: add linebreak to doctests
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 6 years ago

Changed commit from a9b43eb to 228529a

videlec commented 6 years ago
comment:6

Note also that // does the job

sage: L.<q> = LaurentPolynomialRing(QQ)
sage: 1 // q
q^-1
sage: parent(_)
Univariate Laurent Polynomial Ring in q over Rational Field
tscrim commented 6 years ago
comment:7

Replying to @videlec:

Replying to @tscrim:

I disagree with 1/x not being a Laurent polynomial (and I think by extension ~x). This is far too common for people to do and then have to deal with rational functions for something that looks like a Laurent polynomial.

People are familiar with integers and

sage: parent(1 / 1)
Rational Field

Here 1 / 1 looks like an integer. Also people are more familiar with polynomials rather than Laurent polynomials and

sage: R = QQ['x']
sage: parent(R.one() / R.one())
Fraction Field of Univariate Polynomial Ring in x over Rational Field

Still R.one() / R.one() looks like an integer.

The problem with this argument is methods. Elements of ZZ and QQ have nearly identical methods so ducktyping will work between the two. However, this is a very different story for (Laurent) polynomials versus rational functions, in part because the latter is very generic. (I think part of the reason why we use the fraction field of the polynomial ring is so there is less ambiguity with things like 1/x versus x-1.)

IIRC, there were some sage-devel discussions on this. If you do want to make this change, then you need to ask on sage-devel.

I am willing to ask on sage-devel as the current situation is a bug with respect to

The parent of `a / b` should only depend on the parent of `a` and `b`.

I can imagine two alternatives

  1. Allow a / b only when b is invertible which is what you seem to suggest. As a consequence 1 / (1 + q) will raise an error and one would have to explicitely convert to the fraction field to make it work.

That is not what I am proposing. I am proposing keeping the current behavior: a special case for monomials in a Laurent polynomial ring.

  1. Make a / b always return a rational fraction. This is the behavior of all parents in Sage (or at least integers and polynomials). This is also what I proposed in the current branch.

I would also appreciate it if you could put back in the linebreaks on the long outputs.

This I can do. The changes in output look ok? They seem to come from the change in ordering for Laurent polynomials.

Thank you. Yes, I do as they are behaving the same as polynomials, in particular, -x < 0 is now consistent. (IMO, having something like X + (-q^2)*Y is a bit ugly, but that is a separate issue.)

The _new_c also does not mitigate the subclass issue. In fact, it makes it worse if I want to do a Python subclass as I now have no way to override _new_c used in those other methods. I would instead use

cdef type t = type(self)
cdef LaurentPolynomial_univariate p = t.__new__(t)

following sage.groups.perm_gps.permgroup_element.PermutationGroupElement._new_c.

I have to investigate how Cython is optimizing this.

I think this is a serious issue that needs to be addressed:

sage: from sage.rings.polynomial.laurent_polynomial import LaurentPolynomial_univariate
sage: class Foo(Parent):
....:     def __init__(self):
....:         Parent.__init__(self, base=QQ, category=Rings())
....:         self._assign_names('x')
....:     def an_element(self):
....:         return self.element_class(self, [3], -3)
....:     def polynomial_ring(self):
....:         return PolynomialRing(QQ, 'x')
....:     def fraction_field(self):
....:         return self.polynomial_ring().fraction_field()
....:     class Element(LaurentPolynomial_univariate):
....:         pass
....:     
sage: F = Foo()
sage: x = F.an_element(); x
3*x^-3
sage: type(x)
<class '__main__.Foo_with_category.element_class'>
sage: x.inverse_of_unit()
1/3*x^3
sage: type(x.inverse_of_unit())
<type 'sage.rings.polynomial.laurent_polynomial.LaurentPolynomial_univariate'>
tscrim commented 6 years ago

Reviewer: Travis Scrimshaw

videlec commented 6 years ago
comment:9

Replying to @tscrim:

Replying to @videlec:

Replying to @tscrim:

I disagree with 1/x not being a Laurent polynomial (and I think by extension ~x). This is far too common for people to do and then have to deal with rational functions for something that looks like a Laurent polynomial.

People are familiar with integers and

sage: parent(1 / 1)
Rational Field

Here 1 / 1 looks like an integer. Also people are more familiar with polynomials rather than Laurent polynomials and

sage: R = QQ['x']
sage: parent(R.one() / R.one())
Fraction Field of Univariate Polynomial Ring in x over Rational Field

Still R.one() / R.one() looks like an integer.

The problem with this argument is methods. Elements of ZZ and QQ have nearly identical methods so ducktyping will work between the two. However, this is a very different story for (Laurent) polynomials versus rational functions, in part because the latter is very generic. (I think part of the reason why we use the fraction field of the polynomial ring is so there is less ambiguity with things like 1/x versus x-1.)

Ducktyping is dangerous

sage: 3.is_unit()
False
sage: (3/1).is_unit()
True

Moreover, having Laurent polynomials behaving in a different way that any other element in Sage is not helping users. They would better learn that a / b has a parent that only depends on parent(a) and parent(b) and not parent by parent special cases.

IIRC, there were some sage-devel discussions on this. If you do want to make this change, then you need to ask on sage-devel.

I am willing to ask on sage-devel as the current situation is a bug with respect to

The parent of `a / b` should only depend on the parent of `a` and `b`.

I can imagine two alternatives

  1. Allow a / b only when b is invertible which is what you seem to suggest. As a consequence 1 / (1 + q) will raise an error and one would have to explicitely convert to the fraction field to make it work.

That is not what I am proposing. I am proposing keeping the current behavior: a special case for monomials in a Laurent polynomial ring.

This is the worse thing ever. If I write

def my_function(a, b):
    c = a / b
    ... # furter work with c

then depending on whether b is a monomial or not I will end up with different objects for c. It is much safer to have predictible output.

  1. Make a / b always return a rational fraction. This is the behavior of all parents in Sage (or at least integers and polynomials). This is also what I proposed in the current branch.

Discussion at: https://groups.google.com/forum/#!topic/sage-devel/48F1MBo-rEc

jdemeyer commented 6 years ago
comment:10

This is a good argument:

Replying to @videlec:

This is the worse thing ever. If I write

def my_function(a, b):
    c = a / b
    ... # furter work with c

then depending on whether b is a monomial or not I will end up with different objects for c. It is much safer to have predictible output.

jdemeyer commented 6 years ago
comment:11

Replying to @videlec:

Replying to @tscrim:

cdef type t = type(self)
cdef LaurentPolynomial_univariate p = t.__new__(t)

following sage.groups.perm_gps.permgroup_element.PermutationGroupElement._new_c.

I have to investigate how Cython is optimizing this.

Cython deals with this quite well. It optimizes both type(self) and t.__new__(t).

tscrim commented 6 years ago
comment:12

Replying to @jdemeyer:

This is a good argument:

Replying to @videlec:

This is the worse thing ever. If I write

def my_function(a, b):
    c = a / b
    ... # furter work with c

then depending on whether b is a monomial or not I will end up with different objects for c. It is much safer to have predictible output.

However, that should break if b is not a monomial and you are expecting a Laurent polynomial for c. If you know that c generically should be a Laurent polynomial, even when b is not a monomial, then you should be using // anyways. So I don't think this is a good reason.

Addendum - Upon re-reading this, I realized that such code should be using //, even when b is a monomial. However, I think it makes the code harder to debug (see comment:13 below).

tscrim commented 6 years ago
comment:13

Replying to @videlec:

Replying to @tscrim:

The problem with this argument is methods. Elements of ZZ and QQ have nearly identical methods so ducktyping will work between the two. However, this is a very different story for (Laurent) polynomials versus rational functions, in part because the latter is very generic. (I think part of the reason why we use the fraction field of the polynomial ring is so there is less ambiguity with things like 1/x versus x-1.)

Ducktyping is dangerous

sage: 3.is_unit()
False
sage: (3/1).is_unit()
True

That is misleading IMO because it is asking two different questions because of the parent. However, I agree that ducktyping for everything does not work, but for most cases, it is fine (e.g., range(4/2)).

Moreover, having Laurent polynomials behaving in a different way that any other element in Sage is not helping users. They would better learn that a / b has a parent that only depends on parent(a) and parent(b) and not parent by parent special cases.

I agree that this is a good general principle, but having such a hard rule will likely turn users away because of the unexpected behaviors. It also forces people to write code that will have logic bugs that are harder to find and detect than errors that will be raised during computations.

videlec commented 6 years ago
comment:14

merge failures because of #25737