Opened 4 years ago

Closed 9 months ago

#19171 closed defect (fixed)

Add a method `divides` to Polynomial

Reported by: bruno Owned by:
Priority: major Milestone: sage-8.4
Component: commutative algebra Keywords: polynomial, division
Cc: Merged in:
Authors: Bruno Grenet Reviewers: Vincent Delecroix
Report Upstream: N/A Work issues:
Branch: 94c0390 (Commits) Commit: 94c0390bbe4c789678ca2b4f3c730dc97420bb3a
Dependencies: #16649, #25277 Stopgaps:

Description

The generic method divides that can be found in src/sage/structure/element.pyx uses quo_rem (via %) to test if an element divides another one: return (x % self) == 0.

For polynomials, depending on the implementations, the method quo_rem may raise an error if the divisor is not monic (see #16649 for more on this).

This ticket aims at implementing a method divides for the class Polynomial, so that it catches the errors in quo_rem to return False when it is needed.

Example of a problematic behavior:

sage: R.<x> = PolynomialRing(ZZ, implementation="FLINT")
sage: p = 2*x + 1
sage: q = R.random_element(10)
sage: p.divides(q)
False
sage: R.<x> = PolynomialRing(ZZ, implementation="NTL")
sage: p = R(p)
sage: q = R(q)
sage: p.divides(q)
Traceback (most recent call last):
...
ArithmeticError: division not exact in Z[x] (consider coercing to Q[x] first)

Change History (51)

comment:1 Changed 4 years ago by bruno

  • Dependencies set to 16649

comment:2 Changed 4 years ago by bruno

  • Dependencies changed from 16649 to #16649

comment:3 Changed 4 years ago by bruno

  • Branch set to u/bruno/t19171_divides_poly

comment:4 Changed 4 years ago by bruno

  • Commit set to 331c7133714180023a7911c02b6bffd919511a68
  • Status changed from new to needs_review

New commits:

98b7ca1Rebased on 6.9beta1
6d0da24One last zero_element removed + correct one doctest
dbfb518Remove useless colons
331c713#19171: New method divides

comment:5 Changed 3 years ago by git

  • Commit changed from 331c7133714180023a7911c02b6bffd919511a68 to 27040e075be1f5da39bd08bb1d0d9774978bd83b

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

27040e0#19171: New method divides

comment:6 Changed 3 years ago by bruno

  • Component changed from basic arithmetic to commutative algebra
  • Milestone changed from sage-6.9 to sage-7.3
  • Type changed from enhancement to defect

comment:7 Changed 3 years ago by git

  • Commit changed from 27040e075be1f5da39bd08bb1d0d9774978bd83b to 9a189c7bdc29953d6ada184f7b0db503c6da4677

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

9a189c7#19171: New method divides

comment:8 Changed 3 years ago by bruno

I added some noise with two forced pushes (due to some mistake I made), but there is really only one commit! Now the ticket passes the tests (as far as I can tell) and is very ready for review!

comment:9 follow-up: Changed 3 years ago by vdelecroix

You would better use coerce_binop from sage.structure.element rather than a manually handled coercion.

Why is this code specific to polynomials?

comment:10 Changed 3 years ago by vdelecroix

  • Reviewers set to Vincent Delecroix
  • Status changed from needs_review to needs_work

comment:11 Changed 3 years ago by git

  • Commit changed from 9a189c7bdc29953d6ada184f7b0db503c6da4677 to ac43c02c86eb8884bc3d7fa4013f78aef1e34803

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

d03a75e#19171: New method divides
ac43c0219171: Use coerce_binop

comment:12 in reply to: ↑ 9 Changed 3 years ago by bruno

  • Status changed from needs_work to needs_review

Replying to vdelecroix:

You would better use coerce_binop from sage.structure.element rather than a manually handled coercion.

Done.

Why is this code specific to polynomials?

I tried to explain the reasons in the description: For polynomials over (say) ZZ, the euclidean division is not well defined. The problem occurs when the leading coefficient is not invertible. There are several ways to deal with this problem, and different implementations in Sage use different conventions (see #16649 for pointers on this). As a result, testing divisibility by invoking self % p == 0 can lead to ArithmeticError for polynomials: But when this exception is raised, we know that self does not divide p.¹ That's why I decided to add this new method which only differs from the generic one by the fact it catches the exception.

I could have changed the generic method instead to catch an ArithmeticError and return False in such case, but though I am pretty confident that it makes sense for polynomials, I am much less confident concerning other rings.

¹ "We know" is maybe a too strong assumption. In #16649, I did changes that ensure that the assumption holds for all current implementations of polynomials in Sage. (This was not the case before.) Yet, future implementations may break this rule. This may imply that we shouldn't simply rely on quo_rem to test divisibility.

comment:13 follow-up: Changed 3 years ago by vdelecroix

  • Status changed from needs_review to needs_work

Since one is always a unit the if self.is_one(): return True would better be inside the block except NotImplementedError:.

Please add a more exotic example like GF(4)['x']['y'].

comment:14 Changed 3 years ago by vdelecroix

And why not using quo_rem instead of catching a potential ArithmeticError?

comment:15 in reply to: ↑ 13 ; follow-up: Changed 3 years ago by bruno

I am not sure I fully understand your latest comments:

Since one is always a unit the if self.is_one(): return True would better be inside the block except NotImplementedError:.

Do you mean adding a try: ... except NotImplementedError: construction? Because if I add this test say right below the test if self.is_unit(): return True, an exception will be raised in the cases where is_unit is not implemented but is_one is.

And why not using quo_rem instead of catching a potential ArithmeticError?

Using p % self is really the same as p.quo_rem(self)[1] so I do not see the difference.

comment:16 in reply to: ↑ 15 ; follow-up: Changed 3 years ago by vdelecroix

Replying to bruno:

I am not sure I fully understand your latest comments:

Since one is always a unit the if self.is_one(): return True would better be inside the block except NotImplementedError:.

Do you mean adding a try: ... except NotImplementedError: construction? Because if I add this test say right below the test if self.is_unit(): return True, an exception will be raised in the cases where is_unit is not implemented but is_one is.

I meant

try:
    if self.is_unit(): return True   # units divide everything
except NotImplementedError:
    if self.is_one(): return True    # if is_unit is not implemented

And why not using quo_rem instead of catching a potential ArithmeticError?

Using p % self is really the same as p.quo_rem(self)[1] so I do not see the difference.

Not exactly. NTL makes a distinction... and you should actually use pseudo_quo_rem when it is there:

sage: R.<x> = PolynomialRing(ZZ, implementation="NTL")
sage: p = 2*x + 1
sage: q = R.random_element(10)
sage: p.quo_rem(q)
Traceback (most recent call last):
...
ArithmeticError: division not exact in Z[x] (consider coercing to Q[x] first)
sage: p.pseudo_quo_rem(q)
(0, 2*x + 1)

comment:17 Changed 3 years ago by git

  • Commit changed from ac43c02c86eb8884bc3d7fa4013f78aef1e34803 to 0906dfd9801e14acab59d28999fe1778d3cc018d

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

0906dfd19171: Add more involved example + move a line

comment:18 in reply to: ↑ 16 Changed 3 years ago by bruno

  • Status changed from needs_work to needs_review

Replying to vdelecroix:

...

I meant

...

Done.

And why not using quo_rem instead of catching a potential ArithmeticError?

Using p % self is really the same as p.quo_rem(self)[1] so I do not see the difference.

Not exactly. NTL makes a distinction... and you should actually use pseudo_quo_rem when it is there:

sage: R.<x> = PolynomialRing(ZZ, implementation="NTL")
sage: p = 2*x + 1
sage: q = R.random_element(10)
sage: p.quo_rem(q)
Traceback (most recent call last):
...
ArithmeticError: division not exact in Z[x] (consider coercing to Q[x] first)
sage: p.pseudo_quo_rem(q)
(0, 2*x + 1)

Not done: I do not think pseudo_quo_rem is applicable, because of the following:

sage: R.<x> = ZZ[] # works the same with NTL
sage: p = x^2 - 1
sage: q = 2*x + 2
sage: p.pseudo_quo_rem(q)
(2*x - 2, 0)

So, the remainder in p.pseudo_quo_rem(q) can be zero even if q does not divide p.

P.S.: I won't be able to work on this ticket until mid-august.

comment:19 Changed 2 years ago by saraedum

  • Status changed from needs_review to needs_work
  • Work issues set to does not apply

comment:20 Changed 2 years ago by git

  • Commit changed from 0906dfd9801e14acab59d28999fe1778d3cc018d to f85a7e06a7591101005f5e61fff674033344e103

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

f85a7e019171: conflict resolution

comment:21 Changed 2 years ago by bruno

  • Status changed from needs_work to needs_review

comment:22 Changed 2 years ago by tscrim

  • Milestone changed from sage-7.3 to sage-8.0
  • Status changed from needs_review to needs_work

Needs rebasing.

comment:23 Changed 22 months ago by git

  • Commit changed from f85a7e06a7591101005f5e61fff674033344e103 to a0a5badfad9472555eb9b66776ae99ba99416047

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

04b0032Merge branch 'develop' into 19171_divides_poly
a0a5badMerge branch 'u/bruno/t19171_divides_poly' of git://trac.sagemath.org/sage into 19171_divides_poly

comment:24 Changed 22 months ago by bruno

  • Status changed from needs_work to needs_review

comment:25 Changed 22 months ago by mmezzarobba

  • Work issues does not apply deleted

Hi,

This looks reasonable to me; I would be willing to set the ticket to positive review. Vincent, do you agree?

comment:26 Changed 22 months ago by vdelecroix

why return (p % self) == 0 instead of return (p % self).is_zero()?

Last edited 22 months ago by vdelecroix (previous) (diff)

comment:27 follow-up: Changed 22 months ago by vdelecroix

How can you be sure that the following

        except ArithmeticError:
            return False                     # if division is not exact

only catches non exact divisions?

comment:28 Changed 22 months ago by vdelecroix

Could you add check on more exotic rings like Zmod(6)['x'], ZZ['x']['y'], GF(2)['x,y']['z']?

comment:29 in reply to: ↑ 27 Changed 22 months ago by vdelecroix

Replying to vdelecroix:

How can you be sure that the following

        except ArithmeticError:
            return False                     # if division is not exact

only catches non exact divisions?

BTW, each except will cost a non trivial amount of time. Isn't there a fastest way to see if the division is exact? (eg calling quo_rem)

comment:30 Changed 22 months ago by vdelecroix

  • Status changed from needs_review to needs_info

comment:31 Changed 22 months ago by mmezzarobba

FWIW, I agree that using is_zero() would be better, but I don't think it matters enough for blocking the ticket. Considering that an ArithmeticError implies that there is no divisibility looks very reasonable to me, though testing divisibility of the leading coefficients would likely be faster.

comment:32 Changed 18 months ago by git

  • Commit changed from a0a5badfad9472555eb9b66776ae99ba99416047 to 5bac9845a140afc3a913ac2afbeac9346387aabd

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

f350d87Merge branch 'u/bruno/t19171_divides_poly' of git://trac.sagemath.org/sage into t/19171/divides
bcc453419171: Add leading coeff test + change '== 0' to is_zero
5bac98419171: Add tests

comment:33 Changed 18 months ago by bruno

  • Status changed from needs_info to needs_work

I've pushed my branch but some doctest fails. Actually I need some advice.

  • The try... except ArithmeticError is needed for PolynomialRing(ZZ, 'x', implementation = 'NTL') since in this implementation quo_rem raises an exception for non-divisible polynomials, even when the leading coefficients divide each other. See 16649 for details.
  • Contrary to what my (broken) commit message might indicate, I did not add a leading coefficient test that would fasten the computation. The problem comes from the Zmod(6)['x'] example: In Zmod(6), computing the remainder raise exceptions. (Test for instance Zmod(6)(1) % Zmod(6)(4).)
  • As mentioned before, there is still a failing doctest, also due to Zmod(6)['x']. We face the same kind problem as before with the 'NTL' implementation of Z[x]: Computing quo_rem in Zmod(6)['x'] is only possible if the leading coefficient of the divisor is a unit (an exception is raised otherwise). The exception is not caught since contrary to the case of ZZ['x'], the raised exception is a ValueError.

I am not sure to see what the best solutions are to make it work. Some possibilities:

  • We make this ticket depend on another one that unifies all the quo_rem implementations for polynomials over non-field, so that we can base this function on the common behavior of these methods.¹
  • We make this ticket depend on another one with more limited goals, such as unifying the errors to ArithmeticError.
  • We add some specific code, to catch the ValueError.

What do you think?

The first solution seem to me unreachable given the endless discussions about the behavior of quo_rem (links in 16649). The second one is maybe more doable. The third one is probably the faster to close this ticket but I am not sure whether there is a good way to implement it.

¹ More generally, there should exist a ticket for unifying the diverse behaviors of polynomial rings, but 1. it is clearly off-topic this ticket and 2. given the difficulty with unifying "content" in another ticket, I am not at all ready to go for this!

comment:34 Changed 13 months ago by vdelecroix

An ArithmeticError in % does not imply not divisible as shown with

sage: R.<x> = Zmod(6)[]
sage: S.<z> = R[]
sage: p = x*z + 1
sage: q = z + 1
sage: (p*q) % p
Traceback (most recent call last):
...
ArithmeticError: Division non exact (consider coercing to polynomials over the fraction field)

The reason being that

sage: Z6 = Zmod(6)
sage: a = Z6(4)
sage: b = Z6(2)
sage: b.divides(a)
True
sage: a / b
Traceback (most recent call last):
...
ZeroDivisionError: Inverse does not exist.

From there three possibilities:

  • "fix" divisions in integer mod ring Zmod so that a / b above does work. Though in this situation there are several possible answers as 4 = 2 x 5 = 2 x 2 mod 6.
  • restrict divides on polynomial whose base rings are integral domains (so that division is well defined: quotient is unique when it exists).

comment:35 Changed 13 months ago by vdelecroix

A fix is provided in #25277 for Zmod ring.

comment:36 Changed 13 months ago by vdelecroix

  • Dependencies changed from #16649 to #16649, #25277
  • Milestone changed from sage-8.0 to sage-8.3

comment:37 follow-up: Changed 12 months ago by bruno

The difficulties with rings such as Zmod(6) let me think that we should indeed restrict to polynomials over integral domains. The problem being that I actually do not know for instance how to test divisibility over Zmod(6).

If we decide this, the problem with PolynomialRing(ZZ, 'x', implementation="NTL") remains since quo_rem cannot be used there without catching an exception, which is not very appealing... Should we convert first the polynomial to polynomials over the fraction field of the base ring?

comment:38 in reply to: ↑ 37 ; follow-up: Changed 11 months ago by vdelecroix

Replying to bruno:

The difficulties with rings such as Zmod(6) let me think that we should indeed restrict to polynomials over integral domains. The problem being that I actually do not know for instance how to test divisibility over Zmod(6).

There is #25277 for that (it is merged in beta3 already). So you can test for division... but you can not divide yet.

If we decide this, the problem with PolynomialRing(ZZ, 'x', implementation="NTL") remains since quo_rem cannot be used there without catching an exception, which is not very appealing... Should we convert first the polynomial to polynomials over the fraction field of the base ring?

It is fine to catch an exception, but it should be clear what it means.

comment:39 in reply to: ↑ 38 ; follow-up: Changed 11 months ago by bruno

Replying to vdelecroix:

Replying to bruno:

The difficulties with rings such as Zmod(6) let me think that we should indeed restrict to polynomials over integral domains. The problem being that I actually do not know for instance how to test divisibility over Zmod(6).

There is #25277 for that (it is merged in beta3 already). So you can test for division... but you can not divide yet.

I meant testing divisibility of polynomials over Zmod(6). I am not sure that it can be reduced to simply testing divisibility of elements of Zmod(6), I think you need to perform some divisions.

comment:40 in reply to: ↑ 39 Changed 11 months ago by vdelecroix

Replying to bruno:

Replying to vdelecroix:

Replying to bruno:

The difficulties with rings such as Zmod(6) let me think that we should indeed restrict to polynomials over integral domains. The problem being that I actually do not know for instance how to test divisibility over Zmod(6).

There is #25277 for that (it is merged in beta3 already). So you can test for division... but you can not divide yet.

I meant testing divisibility of polynomials over Zmod(6). I am not sure that it can be reduced to simply testing divisibility of elements of Zmod(6), I think you need to perform some divisions.

A far from efficient method consist in

  1. testing divisibility of leading coefficients b_m | a_n
  2. if yes, test all possible quotients q so that a_n = q b_m

But for that, you need a method that return all quotients. Hence my remark.

comment:41 follow-up: Changed 11 months ago by bruno

So do you agree with me that we should better abandon non-integral rings (at least for now)?

comment:42 in reply to: ↑ 41 Changed 11 months ago by vdelecroix

Replying to bruno:

So do you agree with me that we should better abandon non-integral rings (at least for now)?

yes

comment:43 Changed 11 months ago by git

  • Commit changed from 5bac9845a140afc3a913ac2afbeac9346387aabd to 213ccd88e2927bf0df6c06ac81dec003af8e28af

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

fbc129cMerge branch 'u/bruno/t19171_divides_poly' of git://trac.sagemath.org/sage into ticket/19171
213ccd8Add degree test and lc test

comment:44 Changed 11 months ago by git

  • Commit changed from 213ccd88e2927bf0df6c06ac81dec003af8e28af to 235735fd03d71cafe397130b16c830c7c5c65314

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

235735f19171: only for integral domains

comment:45 Changed 11 months ago by bruno

Is the proposed solution ok for you?

comment:46 Changed 11 months ago by vdelecroix

You are not doing the right thing on non-integral domains. You should test is_integral_domain much earlier. With the ticket applied

sage: R.<x> = Zmod(6)[]
sage: (2*x +1).divides(3)
False
sage: (2*x + 1) * 3 == 3
True

(testing divisibility of leading coefficient is not even valid... I was wrong)

You need a line break after TESTS::

It should be documented that this method only works for integral domains. And the examples with Zmod(6) (raising error) should be in the EXAMPLES block rather than TESTS.

comment:47 Changed 11 months ago by git

  • Commit changed from 235735fd03d71cafe397130b16c830c7c5c65314 to 94c0390bbe4c789678ca2b4f3c730dc97420bb3a

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

94c039019171: Divides implemented only for integral domains

comment:48 Changed 11 months ago by bruno

  • Status changed from needs_work to needs_review

I've followed your advice, and changed my mind: My idea was to give an answer (for non-integral domains) as long as I am able to do so. The problem with the divisibility of leading coefficients made me adopt the opposite viewpoint: Divisibility is not implemented for non-integral domains, that's it. I am not sure that an implementation which only works in some uninteresting cases would be helpful to anyone!

comment:49 Changed 10 months ago by vdelecroix

  • Milestone changed from sage-8.3 to sage-8.4

update milestone 8.3 -> 8.4

comment:50 Changed 9 months ago by vdelecroix

  • Status changed from needs_review to positive_review

ok

comment:51 Changed 9 months ago by vbraun

  • Branch changed from u/bruno/t19171_divides_poly to 94c0390bbe4c789678ca2b4f3c730dc97420bb3a
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.