Opened 5 years ago

Closed 5 years ago

#17671 closed defect (fixed)

gcd and xgcd over fields, PID and UFD

Reported by: vdelecroix Owned by:
Priority: major Milestone: sage-6.5
Component: basic arithmetic Keywords:
Cc: Merged in:
Authors: Vincent Delecroix Reviewers: Bruno Grenet
Report Upstream: N/A Work issues:
Branch: 2e0d79e (Commits) Commit: 2e0d79e1c767b471dc143ab55cd5490de342bae0
Dependencies: #17675 Stopgaps:

Description (last modified by vdelecroix)

As reported on sage-support, gcd and xgcd disagree on fraction fields

sage: gcd(6/1,2/1)
2
sage: xgcd(6/1,2/1)
(1, 1/6, 0)

In this ticket we:

  1. Ensure that gcd and xgcd are compatible over quotient fields by implementing the xgcd in the appropriate place. In particular, with that branch we have
    sage: xgcd(6/1,2/1)
    (2, 0, 1)
    sage: xgcd(5/2, 3/4)
    (1/4, -1/2, 2)
    
  2. We fix the xgcd for trivial cases in rings.polynomial.polynomial_element_generic.Polynomial_generic_field in order that gcd and xgcd agree for them.
  3. We introduce a generic test _test_gcd_vs_xgcd in the category PrincipalIdealDomains to ensure the compatibility of gcd and xgcd.
  4. Modify the previous behavior for real numbers
    sage: gcd(2.0, 4.0)
    2
    sage: gcd(2.0, 4.0).parent()
    Integer Ring
    
    The only difference is that now, Sage returns floating point with the same parent
    sage: gcd(2.0, 4.0)
    2.00000000000000
    sage: gcd(2.0, 4.0).parent()
    Real Field with 53 bits of precision
    
    And the same kind of behavior will also hold for lcm and xgcd.

Note:

There is a method xgcd implemented for univariate polynomial over ZZ. This is not a proper name since it has nothing to do with gcd This might not be the best name since the first term is not the gcd in general (see this sage-devel thread and #17674)

Change History (33)

comment:1 Changed 5 years ago by vdelecroix

I am annoyed with the current implementation... it is written in the generic gcd of categories.fields that

For fields of characteristic zero (i.e., containing the
integers as a sub-ring), evaluation in the integer ring is
attempted. This is for backwards compatibility.::
    sage: gcd(6.0,8); gcd(6.0,8).parent()
    2
    Integer Ring

This is completely crazy! And moreover different from what is in gcd!!

comment:2 Changed 5 years ago by vdelecroix

  • Dependencies set to #17673

comment:3 Changed 5 years ago by vdelecroix

  • Branch set to u/vdelecroix/17671
  • Commit set to d8a03d7fa1dd0e401adbae8e2a16ccbfb1902392
  • Status changed from new to needs_review

New commits:

7d2009ctrac #17673: __module__/__name__ for coerce_binop
16d942dtrac #17671: compatibility of gcd and xgcd
d8a03d7trac #17671: fix doctests

comment:4 Changed 5 years ago by vdelecroix

  • Description modified (diff)

comment:5 Changed 5 years ago by vdelecroix

  • Summary changed from Implement gcd for fraction fields to Implement xgcd for fraction fields

comment:6 follow-up: Changed 5 years ago by bruno

In src/sage/categories/quotient_fields.py, lines 262-270:

            except (AttributeError, NotImplementedError, TypeError, ValueError):
                zero = self.parent.zero()
                one  = self.parent.one()
                if self == zero:
                    return (other, zero, one)
                elif other == zero:
                    return (self, one, zero)
                else:
                    return (zero, zero, zero)

I would have returned one as soon as one of self and other is nonzero, and zero only when none is nonzero. The current code is not consistent (if I am not mistaken) with the implementation of gcd (lines 125-128).

comment:7 in reply to: ↑ 6 Changed 5 years ago by vdelecroix

Replying to bruno:

In src/sage/categories/quotient_fields.py, lines 262-270: ... I would have returned one as soon as one of self and other is nonzero, and zero only when none is nonzero.

You are totally right! Thanks!

Vincent

comment:8 Changed 5 years ago by git

  • Commit changed from d8a03d7fa1dd0e401adbae8e2a16ccbfb1902392 to 401193ce270c603dc80ebc1028ff46d0eb662c6c

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

f22ac97trac #17671: (review) fix gcd/xgcd in QuotientFields
401193ctrac #17671: fix xgcd of polynomial_field_generic

comment:9 follow-up: Changed 5 years ago by bruno

  • In src/sage/rings/arith.py, lines 1901-1910, the introducing sentence should to my mind be: Here is an example of a xgcd for two polynomials over the integers, where the linear combination is not the gcd but the gcd multiplied by the resultant::.
        Here is an example of a xgcd for two polynomials over the integers, where the linear
        combination is not the gcd but the resultant::
    
            sage: R.<x> = ZZ[]
            sage: gcd(2*x*(x-1), x^2)
            x
            sage: xgcd(2*x*(x-1), x^2)
            (2*x, -1, 2)
            sage: (2*(x-1)).resultant(x)
            2
    
  • I do not understand exactly what you want to do with your function _test_gcd_vs_xgcd(): As discussed, it is not clear that the results should always be the same.
  • To answer your questions:
    1. As I've said on sage-devel, I am clearly in favor of gcd(RR(6), RR(3)).parent() to return Real Field with 53 bits of precision rather than Integer Ring.
    2. I would keep the name xgcd though, with a clear documentation.

New commits:

f22ac97trac #17671: (review) fix gcd/xgcd in QuotientFields
401193ctrac #17671: fix xgcd of polynomial_field_generic

comment:10 in reply to: ↑ 9 ; follow-up: Changed 5 years ago by vdelecroix

Replying to bruno:

  • In src/sage/rings/arith.py, lines 1901-1910, the introducing sentence should to my mind be:

right

  • I do not understand exactly what you want to do with your function _test_gcd_vs_xgcd(): As discussed, it is not clear that the results should always be the same.

The test is only for principal ideal domains. Up to now I encountered two cases where the convention are a bit special:

  • quotient fields: the gcd/xgcd is compatible with the ring (managed by the category)
  • polynomial over fields: the gcd/xgcd always return monic polynomials (managed by the methods in rings.polynomial.polynomial_element_generic)

If there is an example of a PID where gcd and xgcd must be different I would be happy to remove my test.

Vincent

comment:11 Changed 5 years ago by git

  • Commit changed from 401193ce270c603dc80ebc1028ff46d0eb662c6c to 1d84c4a9160533338cb19a8df2351cd0fc9d79d0

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

1d84c4atrac #17671: (review) clearer doc in rings.arith.xgcd

comment:12 Changed 5 years ago by vdelecroix

  • Description modified (diff)
  • Summary changed from Implement xgcd for fraction fields to gcd and xgcd over fields, PID and UFD

comment:13 in reply to: ↑ 10 Changed 5 years ago by bruno

Replying to vdelecroix:

  • I do not understand exactly what you want to do with your function _test_gcd_vs_xgcd(): As discussed, it is not clear that the results should always be the same.

The test is only for principal ideal domains.

Oh, right! Then I agree with your test.

Everything you've committed so far looks good to me (and passes the tests). There are two things left:

  1. The parent of gcd(3.0,6.0): This ticket seems the right one to change this.
  2. The name of xgcd for polynomials over UFD: There does not seem to be a consensus yet, so I would be in favor of letting this question to a later ticket so that this ticket can get a positive review soon.

comment:14 Changed 5 years ago by vdelecroix

Nope. There seems to be a regression with respect to timings (at least I do not understand what is happening). If I run the doctest on the file coercion_and_categories.rst in $SAGE_ROOT/src/doc/en/thematic_tutorials/ then I end up with an incredibly long TestSuite(P).run(). If I instead replace it with TestSuite(P).run(skip='_test_gcd_vs_xgcd') everything looks fine.

Moreover, in sage.rings.function_field.function_field.py, line 62, the test

TestSuite(M).run(skip = '_test_gcd_vs_xgcd')  # long time (52s on sage.math, 2012)

is much longer for me. And if I remove the skip then it is even longer!

This is a critical issue.

comment:15 Changed 5 years ago by vdelecroix

  • Dependencies changed from #17673 to #17673, #17675

I finally find out the problem in the tests:

  • in function_field.py the timings are similar to what is here in sage-6.5.beta6. This is bad but at least not because of me
  • in coercions_and_categories.rst the problem comes from xgcd over ZZ['x'] (see #17675)

I will add a small commit to fix some doctests.

comment:16 Changed 5 years ago by git

  • Commit changed from 1d84c4a9160533338cb19a8df2351cd0fc9d79d0 to af1f7c61c4dbcf5aa79e95e1c7b62f75b0e18f3f

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

1f84c01trac #17671: revert changes in function_field.py
af1f7c6trac #17671: fix doctest in coercion_and_categories.rst

comment:17 Changed 5 years ago by git

  • Commit changed from af1f7c61c4dbcf5aa79e95e1c7b62f75b0e18f3f to 3275f3dd6e3322b6a7e14f006fa832ebbffd622f

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

cc32b2btrac #17671: gcd(2.0, 6.0) is now 1.0
d05ac70trac #17675: fix xgcd of constant integer polynomials
f0d7ca3trac #17675: doc typo
cf9d62b17675: reviewer's patch: fix typo
089b629merge #17675
91222b2trac #16761: force the parent of gcd/lcm/xgcd
3275f3dtrac #16761: fix typo in doc

comment:18 Changed 5 years ago by git

  • Commit changed from 3275f3dd6e3322b6a7e14f006fa832ebbffd622f to c6db5f832b23d29026019de5306472a25e591279

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

02462bbtrac #17671: force the parent of gcd/lcm/xgcd
c6db5f8trac #17671: fix typo in doc

comment:19 Changed 5 years ago by vdelecroix

Salut Bruno,

Sorry I had a typo in the two last commit messages so I had to do a forced push.

Now everything looks all right. Note that before setting this one to positive review, we need to wait for the trivial #17673 to be in positive review (!hint!). Thanks to Ralf #17675 is already in positive review.

In the last commits:

  • I changed the behavior for real numbers as you and I wanted (nobody answered on the sage-devel thread).
  • I made one change to _test_gcd_vs_xgcd in order to check whether the parent of gcd is really what we expect and not an element of another ring.

Vincent

comment:20 Changed 5 years ago by vdelecroix

  • Description modified (diff)

comment:21 follow-up: Changed 5 years ago by bruno

Salut Vincent,

I looked at #17673, but I am unable to review it since it goes beyond my knowledge in Sage... Or at least it will take some time for me to understand it and then review.

Concerning the real number, I would have thought of an implementation acting as follows:

sage: gcd(2.0,6.0)
2.00000000000000

I do not know if this makes sense.

comment:22 in reply to: ↑ 21 ; follow-up: Changed 5 years ago by vdelecroix

Replying to bruno:

Salut Vincent,

I looked at #17673, but I am unable to review it since it goes beyond my knowledge in Sage... Or at least it will take some time for me to understand it and then review.

Concerning the real number, I would have thought of an implementation acting as follows:

sage: gcd(2.0,6.0)
2.00000000000000

Why? Would you like it to be nice only over ZZ and not over QQ? And this is where it becomes impossible. Within Sage any floating point can be converted to QQ...

sage: QQ(1.1231938951)
155072989/138064309

So the only way I would agree with your proposition is to have also

sage: gcd(RR(5/2), RR(3/2)) == RR(gcd(5/2, 3/2))

which is currently False as well.

Anyway mixing algebraic functions (like gcd) with floating point number is often a bad idea. So at least the proposed implementation tells you that something is wrong.

Vincent

comment:23 in reply to: ↑ 22 Changed 5 years ago by bruno

Replying to vdelecroix:

Why? Would you like it to be nice only over ZZ and not over QQ? And this is where it becomes impossible.

First of all, your current solution is mathematically perfectly correct and I am quite fine with it. This said, let me expose the reasons of my suggestion.

With my suggestion, "integral floating point numbers" are treated differently than other floating point numbers. Reasons:

  1. The gcd function coincide with the mathematical GCD function for integers only;
  2. Every floating point number can be seen as a rational, while not every floating point number can be seen as an integer.

Also, RR and QQ are treated differently. Reasons:

  1. For inputs in QQ, we know that they are rational numbers, so we deviate from the mathematical definition to give a more meaningful answer and the definition we use "lives" in QQ;
  2. For inputs in RR that are "integral", we make an assumption (the floating point numbers we have at hand represent integers) and then use the standard mathematical definition;
  3. For inputs in RR that are not "integral", not returning 1.0000000 requires to make an assumption AND to deviate from the mathematical definition: that is, we would need to "move" from RR to QQ and then to use a non-standard definition.

Finally, if a user asks for the GCD of two "integral floating point numbers", I have no difficulty to understand what s/he means (since the GCD of the integers they represent is a well-defined mathematical function). It is much more difficult to me to interpret the will of a user that asks for the GCD of two "rational floating point numbers", and this (also) justifies to treat differently the two situations.

Anyway mixing algebraic functions (like gcd) with floating point number is often a bad idea. So at least the proposed implementation tells you that something is wrong.

I agree with that.

comment:24 Changed 5 years ago by vdelecroix

All right. Let us do it your way.

(actually mine would be to not define floating point as a field ;-)

comment:25 Changed 5 years ago by git

  • Commit changed from c6db5f832b23d29026019de5306472a25e591279 to 0f536c92b61f088375410417bc790fee6de8dc10

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

db90fddtrac #17671: (review) fix gcd/xgcd in QuotientFields
673bca9trac #17671: fix xgcd of polynomial_field_generic
dd86cb7trac #17671: (review) clearer doc in rings.arith.xgcd
5bc4805trac #17671: revert changes in function_field.py
994196ftrac #17671: fix doctest in coercion_and_categories.rst
c073eb9trac #17671: gcd(2.0, 6.0) is now 1.0
9714f6atrac #17671: force the parent of gcd/lcm/xgcd
9bbb25ctrac #17671: fix typo in doc
5d9f8d2trac #17671: merge #17675
0f536c9trac #17671: fix float definition of gcd/lcm/xgcd

comment:26 Changed 5 years ago by vdelecroix

  • Dependencies changed from #17673, #17675 to #17675

I removed #17673 from the dependencies... I had to do a forced push to remove the commit coming from that branch.

Vincent

comment:27 Changed 5 years ago by bruno

It seems that I get the same problem you had in #17673, thus the tests in src/sage/rings/polynomial/polynomial_quotient_ring.py don't pass. From line 281, there is a test that looks like the following (I removed irrelevant parts):

sage: P.<x> = QQ[]
sage: Q = P.quotient(x^2+2)
sage: Q in Fields()
True
sage: Q.an_element().gcd._module_
'sage.categories.fields'

but the test doesn't pass, and I get AttributeError: 'sage.structure.element.NamedBinopMethod' object has no attribute '__module__' instead of 'sage.categories.fields'.

comment:28 Changed 5 years ago by vdelecroix

Yep. I forgot to change that. And the modification of gcd/lcm for floating point also causes some trouble for the symbolic ring. A commit is comming.

Vincent

comment:29 Changed 5 years ago by git

  • Commit changed from 0f536c92b61f088375410417bc790fee6de8dc10 to 2e0d79e1c767b471dc143ab55cd5490de342bae0

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

2e0d79etrac #17671: fix two doctests

comment:30 Changed 5 years ago by vdelecroix

  • Description modified (diff)

comment:31 Changed 5 years ago by bruno

  • Reviewers set to Bruno Grenet
  • Status changed from needs_review to positive_review

All tests pass! LGTM.

comment:32 Changed 5 years ago by vdelecroix

Cool! Thanks for the review.

comment:33 Changed 5 years ago by vbraun

  • Branch changed from u/vdelecroix/17671 to 2e0d79e1c767b471dc143ab55cd5490de342bae0
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.