Opened 9 years ago

Closed 8 years ago

#10771 closed enhancement (fixed)

gcd and lcm for fraction fields

Reported by: SimonKing Owned by: AlexGhitza
Priority: major Milestone: sage-4.7.2
Component: basic arithmetic Keywords: gcd lcm fraction fields
Cc: burcin Merged in: sage-4.7.2.alpha0
Authors: Simon King Reviewers: Marco Streng, Mariah Lenox
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Description (last modified by SimonKing)

At sage-devel, the question was raised whether it really is a good idea that the gcd in the rational field should return either 0 or 1.

Since any non-zero element of QQ qualifies as gcd of two non-zero rationals, it should be possible to define gcd and lcm, so that gcd(x,y)*lcm(x,y)==x*y holds for any rational numbers x,y, and so that gcd(QQ(m),QQ(n))==gcd(m,n) and lcm(QQ(m),QQ(n))==lcm(m,n) for any two integers m,n.

Moreover, it should be possible to provide gcd/lcm for any fraction field of a PID: Note that currently gcd raises a type error for elements of Frac(QQ['x']).

The aim is to implement gcd and lcm as ElementMethods of the category QuotientFields().

Approach

Let R be an integral domain, assume that it provides gcd and lcm, and let F be its fraction field. Since R has gcd, we can assume that x.numerator() and x.denominator() are relatively prime for any element x of F.

Then, define

gcd(x,y) = gcd(x.numerator(),y.numerator())/lcm(x.denominator(),y.denominator())
lcm(x,y) = lcm(x.numerator(),y.numerator())/gcd(x.denominator(),y.denominator())

Benefits

If that approach is mathematically sober, we obtain the following equalities up to units in R:

  • gcd(x,y)*lcm(x,y)==x*y, for any x,y in F, provided that the equality holds for any x,y in R.
  • gcd(F(x),F(y))==gcd(x,y) and lcm(F(x),F(y))==lcm(x,y) for any x,y in R.

Attachments (1)

trac10771-lcm_gcd_fraction_fields.patch (19.7 KB) - added by SimonKing 8 years ago.
Implement gcd and lcm for general fields, with the special case of the fraction field of a unique factorization domain

Download all attachments as: .zip

Change History (39)

comment:1 Changed 9 years ago by SimonKing

  • Description modified (diff)

comment:2 Changed 9 years ago by SimonKing

  • Authors set to Simon King
  • Status changed from new to needs_review
  • Type changed from PLEASE CHANGE to defect

Without the patch, we had:

sage: gcd(2/1,4)  # does not restrict to ZZ
1
sage: lcm(2/1,4)  # Bug
Traceback (most recent call last):
...
TypeError: Argument 'other' has incorrect type (expected sage.rings.rational.Rational, got sage.rings.integer.Integer)
sage: R.<x> = QQ[]
sage: lcm(1/(x+1),1/(x+1)^2) # note that the error message names gcd, not lcm!
Traceback (most recent call last):
...
TypeError: unable to find gcd of 1/(x + 1) and 1/(x^2 + 2*x + 1)
sage: gcd(1/(x+1),1/(x+1)^2)
Traceback (most recent call last):
...
TypeError: unable to find gcd of 1/(x + 1) and 1/(x^2 + 2*x + 1)
sage: gcd(int(2),2/1)
2
sage: gcd(2,2/1) # gcd of ints and integers are different
1

With the patch, one has

sage: gcd(2/1,4)
2
sage: lcm(2/1,4)
4
sage: R.<x> = QQ[]
sage: lcm(1/(x+1),1/(x+1)^2)
1/(x + 1)
sage: gcd(1/(x+1),1/(x+1)^2)
1/(x^2 + 2*x + 1)
sage: gcd(int(2),2/1)
2
sage: gcd(2,2/1)
2

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

  • Reviewers set to Marco Streng

A couple of comments:

  • I don't see the point of keeping a special function gcd_rational(self, other, **kwds): that returns a GCD in the set {0,1}. Why only for QQ? Why at all? Also, the difference between gcd and gcd_rational is not explained by the word "rational". Perhaps gcd_zero_one would be a more informative name.
  • Everything generalizes from principal ideal domains to unique factorization domains (UFD's) (such as multivariate polynomial rings over unique factorization domains) as long as they have gcd and lcm methods implemented. Why not write "unique factorization domain" in the documentation instead of "principal ideal domain"?
  • Suppose I have a fraction field F of a ring of type R that does not have lcm and gcd methods, or these methods exist, but raise other kinds of errors, e.g. because the ring is not a UFD or the methods have not been implemented. Let a and b be elements of F. Then a and b have a gcd in F because F is a field, so I would expect a.gcd(b) to return something (anything basically). After applying your patch, if I do a.gcd(b), it is very confusing to get an AttributeError: 'RElement' object has no attribute 'gcd': I'm not interested in gcd's of RElements, only of (Fraction)FieldElements. You could put your entire gcd and lcm code between try: and except (AttributeError, NotImplementedError, TypeError, ValueError): to return the same 0 or 1 that gcd_rational would (which is a mathematically correct gcd in F after all).
  • trac has a t too much in line 1610 of ring.pyx: the case of the rational field. However, since tract ticket #10771,
  • you write "quotient field" in the documentation. You could write "fraction field" to avoid confusion with quotient rings, which may be fields. This makes it more clear that you refer to the mathematical counterpart of Sage's "FractionField"?
  • All tests passed! with -long, well documented, looks good.

comment:4 Changed 9 years ago by lftabera

Related tickets: #9819, #10459

comment:5 Changed 9 years ago by dsm

Trivial: typo "tract" for "trac" in patch.

comment:6 in reply to: ↑ 3 Changed 9 years ago by SimonKing

Hi Marco,

Replying to mstreng:

A couple of comments:

  • I don't see the point of keeping a special function gcd_rational(self, other, **kwds): that returns a GCD in the set {0,1}. Why only for QQ? Why at all? Also, the difference between gcd and gcd_rational is not explained by the word "rational". Perhaps gcd_zero_one would be a more informative name.

Perhaps I am too 'conservative': I would rather rename something than to delete it.

But I suppose it would be better to implement a gcd and an lcm as element methods of Fields() -- see #9819. In that way, one could easily provide a correct gcd/lcm behaviour for all fields, and for fraction fields of PID one would still obtain a gcd/lcm behaviour that is not only correct but nice.

  • Everything generalizes from principal ideal domains to unique factorization domains (UFD's) (such as multivariate polynomial rings over unique factorization domains) as long as they have gcd and lcm methods implemented. Why not write "unique factorization domain" in the documentation instead of "principal ideal domain"?

Good question. I did it by accident. However, it matches the sad fact that the categories are a little askew here:

sage: PrincipalIdealDomains().is_subcategory(UniqueFactorizationDomains())
False  # shouldn't that be "True"?
sage: PrincipalIdealDomains().is_subcategory(GcdDomains())
True
sage: UniqueFactorizationDomains().is_subcategory(GcdDomains())
True

Perhaps one should write "fraction field of an integral domain with gcd and lcm"? Because that's what is duck typed.

  • Suppose I have a fraction field F of a ring of type R that does not have lcm and gcd methods, or these methods exist, but raise other kinds of errors, e.g. because the ring is not a UFD or the methods have not been implemented. Let a and b be elements of F. Then a and b have a gcd in F because F is a field, so I would expect a.gcd(b) to return something (anything basically). After applying your patch, if I do a.gcd(b), it is very confusing to get an AttributeError: 'RElement' object has no attribute 'gcd': I'm not interested in gcd's of RElements, only of (Fraction)FieldElements. You could put your entire gcd and lcm code between try: and except (AttributeError, NotImplementedError, TypeError, ValueError): to return the same 0 or 1 that gcd_rational would (which is a mathematically correct gcd in F after all).

If one adds a gcd and lcm for field elements (#9819) returning 0 or 1, then your suggestion certainly makes sense.

  • trac has a t too much in line 1610 of ring.pyx: the case of the rational field. However, since tract ticket #10771,
  • you write "quotient field" in the documentation. You could write "fraction field" to avoid confusion with quotient rings, which may be fields. This makes it more clear that you refer to the mathematical counterpart of Sage's "FractionField"?

Thanks! I am about to post a new patch, and I hope that Luis as author of #9819 does not mind if I add the case of arbitrary fields to my patch.

comment:7 follow-up: Changed 9 years ago by SimonKing

I just updated my patch, and I hope it addresses all your remarks. The typos are fixed, the old gcd of the rational field is now removed, and: I also added gcd/lcm for general fields. I am sorry that this makes #9819 a duplicate.

And I hope that I did not confuse gcd and lcm in the following setting:

sage: GF(2)(1).gcd(GF(2)(1))
1
sage: GF(2)(1).gcd(GF(2)(0))
1
sage: GF(2)(0).gcd(GF(2)(0))
0
sage: GF(2)(1).lcm(GF(2)(0))
0
sage: GF(2)(1).lcm(GF(2)(1))
1

That's implemented as element methods for the category of Fields(). Somehow, the category framework is cool, isn't it?

comment:8 follow-up: Changed 9 years ago by SimonKing

Oops, I forgot one aspect: If the gcd/lcm of the base ring raises an error, it should be caught. Implementing it now...

comment:9 in reply to: ↑ 8 Changed 9 years ago by SimonKing

Replying to SimonKing:

Oops, I forgot one aspect: If the gcd/lcm of the base ring raises an error, it should be caught. Implementing it now...

Here it is:

sage: R.<q> = ZZ.extension(x^2+5)
sage: gcd(q,q)
1
sage: gcd(q,0)
1
sage: gcd(R.zero(),0)
0
sage: lcm(q,q)
1
sage: lcm(q,0)
0

Correct me if I am wrong, but I think that now every complaint is addressed.

Best regards,

Simon

comment:10 follow-up: Changed 9 years ago by mstreng

  • Status changed from needs_review to needs_work

New complaints, sorry.

  • I don't really like this:
    sage: M = Matrix(GF(3),2,2,[1,2,0,1]); M
    [1 2]
    [0 1]
    sage: GF(2)(1).lcm(M)
    1
    
    I think you should first try to coerce other to self.parent() in gcd and lcm for Field elements.
  • Line 64 of quotient_fields.py "both GCD and LCM, it is possible to be a bit more specific" uses capital letters GCD and LCM, while the methods that it (quasi-)refers to are lower case: AttributeError: 'sage.rings.rational.Rational' object has no attribute 'GCD'
  • Next line: "define the GCD uniquely up to a unit in the base ring". Write "of" instead of "in" so that it is clear that you are not saying that the GCD is in the base ring.
  • Your tests with R.<q> = ZZ.extension(x^2+5) is unrelated to fraction field: q is 1
    sage: R.<q> = ZZ.extension(x^2+5)
    sage: g = R.gens(); g
    [1, q]
    sage: q
    1
    sage: g[1]
    q
    sage: g[1] == q
    False
    
    All your examples involving q get coerced to ZZ by the global GCD and LCM, so you should find other examples. This seems to work:
    sage: gcd(g[1],1)
    TypeError: unable to find gcd of q and 1
    sage: gcd(g[1]/1,1)
    1
    
    (Correct behaviour both times)

Time to catch a bus, may have more comments later.

comment:11 in reply to: ↑ 7 ; follow-up: Changed 9 years ago by lftabera

Replying to SimonKing:

I just updated my patch, and I hope it addresses all your remarks. The typos are fixed, the old gcd of the rational field is now removed, and: I also added gcd/lcm for general fields. I am sorry that this makes #9819 a duplicate.

Please.

And I hope that I did not confuse gcd and lcm in the following setting:

sage: GF(2)(1).gcd(GF(2)(1))
1
sage: GF(2)(1).gcd(GF(2)(0))
1
sage: GF(2)(0).gcd(GF(2)(0))
0
sage: GF(2)(1).lcm(GF(2)(0))
0
sage: GF(2)(1).lcm(GF(2)(1))
1

That's implemented as element methods for the category of Fields(). Somehow, the category framework is cool, isn't it?

I do not full understand this, but I am happy that works.

Some thoughts:

For QQ and the like, could it be that gcd and lcm should only take care of coercing to a good setting and then the real algorithm should be in _lcm, _gcd? Note that QQ already has ._gcd and ._lcm so these methods has to be taken into account.

For generic Fields, It should appear in the documentation that the methods return 0 if both arguments are zero and a non-zero element otherwise. The user should not suppose that the non-zero element is actually one, since this is changed by subclasses.

comment:12 in reply to: ↑ 11 Changed 9 years ago by SimonKing

Hi Luis!

Replying to lftabera:

That's implemented as element methods for the category of Fields(). Somehow, the category framework is cool, isn't it?

I do not full understand this, but I am happy that works.

Then let me explain:

Any category (in Sage) comes with a parent class and an "element class". That very often is an abstract thing, but you can inherit from it, and in particular you can implement methods for it.

If things are properly done, then a parent P is declared as an object in some category, and its class automatically(!) inherits from the parent class of its category:

sage: QQ.__class__
<class 'sage.rings.rational_field.RationalField_with_category'>
sage: issubclass(QQ.__class__,QuotientFields().parent_class)
True

If things are very properly done, then the parent has an attribute Element, which is supposed to be a class. That class will (again automatically) be translated into a subclass of the element class of the parent's category, and is available as the attribute element_class:

sage: m = matrix(GF(2),[[1]])
sage: G = MatrixGroup([m])
sage: G.Element
<class 'sage.groups.matrix_gps.matrix_group_element.MatrixGroupElement'>
sage: G.__class__
<class 'sage.groups.matrix_gps.matrix_group.MatrixGroup_gens_finite_field_with_category'>
sage: G.element_class
<class 'sage.groups.matrix_gps.matrix_group_element.MatrixGroup_gens_finite_field_with_category.element_class'>
sage: issubclass(G.element_class,G.Element)
True
sage: issubclass(G.element_class,G.category().element_class)
True
sage: isinstance(G.random_element(),G.category().element_class)
True

So, if things are very properly done, then a method defined as an "element method" of a category, will be available by looking up the method resolution order, since the elements of the parent are actual instances of the category's element class.

Sometimes things are done properly, but not very properly:

sage: QQ.category()
Category of quotient fields
sage: isinstance(QQ, QQ.category().parent_class)
True
sage: hasattr(QQ,'Element')
False
sage: isinstance(QQ.one(), QQ.category().element_class)
False

Hence, the methods that I defined for the element class of the "category of quotient fields" is not available to elements of QQ by simply looking up the method resolution order. But: There also is a __getattr__ method implemented for sage.structure.element.Element, and this can access the element methods of the category of the parent of an element (uff!) even if the element is no instance of the "proper" element class.

This is why the new gcd works for the rationals, but it's slower than with the method resolution order, and this is why I also did not remove the existing lcm for rationals.

Some thoughts:

For QQ and the like, could it be that gcd and lcm should only take care of coercing to a good setting and then the real algorithm should be in _lcm, _gcd? Note that QQ already has ._gcd and ._lcm so these methods has to be taken into account.

That's a nice observation. By searching a little, I found that _gcd is used in the method sage.structure.element.PrincipalIdealDomainElement.gcd.

However, there are only very few classes that define a _gcd method: sage.structure.element.EuclideanDomainElement, sage.structure.element.FieldElement, sage.rings.rational.Rational, and then there are hits in ...polynomial_element_generic and ...polynomial_modn_dense_ntl.

In other words, I guess that the _gcd method is an artefact of earlier attempts to reflect mathematics in the class hierarchy -- this should now be done in the category framework. _gcd will only play a role if

  1. the class of an element derives from PrincipalIdealDomainElement, and
  1. the class of an element additionally defines _gcd.

It seems to me that there is precisely one class that directly inherits from PrincipalIdealDomainElement, and there is precisely one class that directly inherits from that class:

sage: search_src("\(PrincipalIdealDomainElement")
structure/element.pxd:83:cdef class EuclideanDomainElement(PrincipalIdealDomainElement):
structure/element.pyx:2370:cdef class EuclideanDomainElement(PrincipalIdealDomainElement):

sage: search_src("\(EuclideanDomainElement")
structure/element.pyx:2362:PY_SET_TP_NEW(EuclideanDomainElement, Element)
rings/integer.pxd:9:cdef class Integer(EuclideanDomainElement):

So, _gcd is only used for the Integer class. I guess it could easily be removed.

For generic Fields, It should appear in the documentation that the methods return 0 if both arguments are zero and a non-zero element otherwise. The user should not suppose that the non-zero element is actually one, since this is changed by subclasses.

Right. I'll take care of that, but now I have a bus to catch...

comment:13 follow-up: Changed 9 years ago by mstreng

  • Did you run a full doctest? I'm getting doctest failures that may have something to do with your patch, as they probably simplify fractions somewhere deep inside in some way. Or my installation is messed up again.
    sage -t -long devel/sage/sage/symbolic/expression.pyx # 1 doctests failed
    sage -t -long devel/sage/sage/symbolic/integration/integral.py # 1 doctests failed
    sage -t -long devel/sage/sage/stats/basic_stats.py # 2 doctests failed
    
  • The documentation of content for multivariate polynomials says "Returns the content of this polynomial. Here, we define content as the gcd of the coefficients in the base ring." Your changed doctest (1 becomes 2) is correct, but perhaps it would be more informative to add f.content().parent(), and to give an example over another field. (Just a suggestion if you are editing the code anyway. Otherwise, don't bother.) More important:
    # Since trac ticket #10771, the gcd in QQ restricts to the 
    # gcd in ZZ. 
    
    I don't think this comment is needed in this part of the code. If you want to include it, it would better fit in the doctest a few lines above it. And I think junk like the following should not be in the code at all: #,integer=self.parent() is ZZ)

comment:14 in reply to: ↑ 13 Changed 9 years ago by SimonKing

Replying to mstreng:

  • Did you run a full doctest? I'm getting doctest failures that may have something to do with your patch, as they probably simplify fractions somewhere deep inside in some way. Or my installation is messed up again.
    sage -t -long devel/sage/sage/symbolic/expression.pyx # 1 doctests failed
    sage -t -long devel/sage/sage/symbolic/integration/integral.py # 1 doctests failed
    sage -t -long devel/sage/sage/stats/basic_stats.py # 2 doctests failed
    

Hum. Apparently I did not run these tests. The problem is that I currently work on various patches, frequently switch between them, and it may happen that I test one of them but forget the other.

Anyway, the failure in "expression.pyx" seems to be related with the patch. I'll try to understand what is happening.

  • The documentation of content for multivariate polynomials says "Returns the content of this polynomial. Here, we define content as the gcd of the coefficients in the base ring." Your changed doctest (1 becomes 2) is correct, but perhaps it would be more informative to add f.content().parent(), and to give an example over another field. (Just a suggestion if you are editing the code anyway. Otherwise, don't bother.)

OK.

More important:

# Since trac ticket #10771, the gcd in QQ restricts to the 
# gcd in ZZ. 

I don't think this comment is needed in this part of the code. If you want to include it, it would better fit in the doctest a few lines above it.

I thought of it as a comment to people who know the original code. It could be that it is a method that uses gcd only internally, so that commenting in the code seems better than in the documentation.

And I think junk like the following should not be in the code at all: #,integer=self.parent() is ZZ)

OK, that was an artefact of editing.

But I think most urgent are the doctest failures. It seems related to gcd for fields that are not fraction fields (such as RR). That would also explain why all text passed with the original version of my patch, whereas the new version fails.

Cheers, Simon

comment:15 Changed 9 years ago by SimonKing

Yep, I understand where the first error comes from.

The original gcd codes says:

    if b is not None:
        if hasattr(a, "gcd"):
            return a.gcd(b, **kwargs)
        else:
            try:
                return ZZ(a).gcd(ZZ(b))
            except TypeError:
                raise TypeError, "unable to find gcd of %s and %s"%(a,b)

Moreover, elements of RR did not have a gcd method. Hence, it was tried to coerce them into ZZ, and so we had

sage: gcd(3.0,6.0)
3

But with my patch, all field elements (including elements of RR) have a gcd method, so, a conversion to ZZ is not attempted.

I am not sure yet how I would solve it: Change the code of `sage.symbolic.expression"? Change the gcd implemented for fields that are no fraction fields, such that a special case is made for inexact fields?

comment:16 Changed 9 years ago by SimonKing

I suggest to do the following, which would generalize the current attempt to convert stuff into ZZ:

            For inexact fields, evaluation in the prime subfield is attempted::

                sage: lcm(15.2,12.0); lcm(15.2,12.0).parent()
                228
                Rational Field
                sage: RR(76/5)
                15.2000000000000
                sage: lcm(76/5,12)
                228

            If this fails, we resort to the default we see above::

                sage: lcm(6.0*CC.0,8*CC.0); lcm(6.0*CC.0,8*CC.0).parent()
                1.00000000000000
                Complex Field with 53 bits of precision

Similarly for gcd.

In that way, the first of the three doctest failures vanishes. But the other two remain problems.

comment:17 follow-up: Changed 9 years ago by SimonKing

The second failure can be triggered as follows:

sage: F(x) = 1/sqrt(2*pi*1^2)*exp(-1/(2*1^2)*(x-0)^2)
sage: G(x) = 1/sqrt(2*pi*n(1)^2)*exp(-1/(2*n(1)^2)*(x-n(0))^2)
sage: ((F(x)-G(x))^2).expand()
/mnt/local/king/SAGE/sage-4.6.2.alpha4/local/bin/sage-sage: Zeile 300: 18366 Speicherzugriffsfehler  sage-ipython "$@" -i

Segmentation fault while expanding. I could imagine that gcd==1 results in an infinite loop.

comment:18 in reply to: ↑ 17 Changed 9 years ago by SimonKing

Replying to SimonKing:

Segmentation fault while expanding. I could imagine that gcd==1 results in an infinite loop.

Ouch, it is worse: It already occurs before the expansion!

comment:19 Changed 9 years ago by SimonKing

The third error boils down to:

Without patch

sage: x = -sqrt(2)-1/5*I
sage: x*x
1/25*(5*sqrt(2) + I)^2

With patch

sage: x = -sqrt(2)-1/5*I
sage: x*x
1/25*(-5*sqrt(2) - I)^2

The segfault in the second doctest problem also comes from multiplication of symbolic expressions. And the example above indicates that it is related with pulling common factors out of a list of expressions.

comment:20 Changed 9 years ago by kcrisman

Just FYI, #8111 is possibly also related.

comment:21 follow-up: Changed 9 years ago by SimonKing

  • Cc burcin added
  • Status changed from needs_work to needs_review

With the new patch, all tests pass for me. The behaviour is now (mutatis mutandis for lcm):

  • If ZZ is a subring of a field F that is not a fraction field then a.gcd(b) first tests if a and b can be converted into ZZ. If it is the case, their gcd in ZZ as element of ZZ is returned. This is not nice, but it is compatible with the "unpatched" behaviour of Sage and is assumed in some other parts of Sage:
    # unpatched!
    sage: gcd(RR(1),RR(1))
    1
    sage: _.parent()
    Integer Ring
    
  • If ZZ is no subring of F or conversion fails then a.gcd(b) returns F(0) or F(1). This is, e.g., the case for gcd(CC(2),CC(2)), since CC(2) can not be converted into ZZ.
  • a.gcd(b) only raises an error if b can not be converted into a.parent().
  • gcd and lcm in fraction fields of a unique factorization domain R restrict to R, and moreover they satisfy gcd(a,b)*lcm(a,b)==a*b up to multiplication with a unit of R.

There is a price to pay, and that's why I put Burcin on Cc:

The distribution of minus signs in symbolic expression changes, apparently since now gcd raises no errors in examples like the following:

sage: (I - 1/3*sqrt(2))^2
1/9*(-sqrt(2) + 3*I)^2
# Without patch, we would get
# 1/9*(sqrt(2) - 3*I)^2

The question is whether a changed minus sign distribution is such a serious thing.

comment:22 in reply to: ↑ 10 Changed 9 years ago by SimonKing

Hi Marco,

sorry, I forgot to explain how I addressed your complaints:

Replying to mstreng:

New complaints, sorry.

  • I don't really like this:
    sage: M = Matrix(GF(3),2,2,[1,2,0,1]); M
    [1 2]
    [0 1]
    sage: GF(2)(1).lcm(M)
    1
    
    I think you should first try to coerce other to self.parent() in gcd and lcm for Field elements.

Yep, that was a mistake. We now have:

sage: M = Matrix(GF(3),2,2,[1,2,0,1]); M
[1 2]
[0 1]
sage: GF(2)(1).lcm(M)
---------------------------------------------------------------------------
ArithmeticError                           Traceback (most recent call last)
...
ArithmeticError: The second argument can not be interpreted in the parent of the first argument. Can't compute the lcm
sage: GF(2)(1).gcd(M)
ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (4, 0))

---------------------------------------------------------------------------
ArithmeticError                           Traceback (most recent call last)

...
ArithmeticError: The second argument can not be interpreted in the parent of the first argument. Can't compute the gcd

  • Line 64 of quotient_fields.py "both GCD and LCM, it is possible to be a bit more specific" uses capital letters GCD and LCM, while the methods that it (quasi-)refers to are lower case:

I changed the doc string to lower case.

  • Next line: "define the GCD uniquely up to a unit in the base ring". Write "of" instead of "in" so that it is clear that you are not saying that the GCD is in the base ring.

Done.

  • Your tests with R.<q> = ZZ.extension(x^2+5) is unrelated to fraction field:

You are right. My intention was to start with R, show that there is no gcd in R, and then expose what happens in Frac(R). The test is changed accordingly.

Best regards,

Simon

comment:23 in reply to: ↑ 21 ; follow-up: Changed 9 years ago by kcrisman

There is a price to pay, and that's why I put Burcin on Cc:

The distribution of minus signs in symbolic expression changes, apparently since now gcd raises no errors in examples like the following:

sage: (I - 1/3*sqrt(2))^2
1/9*(-sqrt(2) + 3*I)^2
# Without patch, we would get
# 1/9*(sqrt(2) - 3*I)^2

The question is whether a changed minus sign distribution is such a serious thing.

What happens with the already existing doctests? That would be what I would be most concerned about. Also, is the infinite loop now gone? IIRC you had fixed that somehow.

Really, in some sense the new version is 'better' because it keeps the minus where the user put it. Do (-I - 1/3*sqrt(2))^2 and (I + 1/3*sqrt(2))^2 do anything different from before?

comment:24 in reply to: ↑ 23 ; follow-up: Changed 9 years ago by SimonKing

Replying to kcrisman:

What happens with the already existing doctests? That would be what I would be most concerned about.

There are two tests (really in fact only one) in sage.stats where I had to change the expected output. The old and new output are equivalent:

sage: std([I, sqrt(2), 3/5]) # result with patch
sqrt(1/450*(-5*sqrt(2) + 10*I - 3)^2 + 1/450*(-5*sqrt(2) - 5*I + 6)^2 + 1/450*(10*sqrt(2) - 5*I - 3)^2)
# Without patch, we had sqrt(1/450*(5*sqrt(2) + 5*I - 6)^2 + 1/450*(5*sqrt(2) - 10*I + 3)^2 + 1/450*(10*sqrt(2) - 5*I - 3)^2)
# Testing equality:
sage: bool(sqrt(1/450*(5*sqrt(2) + 5*I - 6)^2 + 1/450*(5*sqrt(2) - 10*I + 3)^2 + 1/450*(10*sqrt(2) - 5*I - 3)^2) == std([I, sqrt(2), 3/5]))
True

There is essentially the same happening with variance([I, sqrt(2), 3/5]).

Also, is the infinite loop now gone? IIRC you had fixed that somehow.

Yep. The problem was that in an old version of my patch I would not try conversion to ZZ but conversion to the prime subfield. A real number with finitely many digits can of course always be converted into RR.prime_subfield(), which is QQ. So, I guess that the segfault came from rationals with very large numerator and denominator. Here is the proof that the segfault has gone:

sage: F(x) = 1/sqrt(2*pi*1^2)*exp(-1/(2*1^2)*(x-0)^2) 
sage: G(x) = 1/sqrt(2*pi*n(1)^2)*exp(-1/(2*n(1)^2)*(x-n(0))^2) 
sage: (F-G)**2 
x |--> 1/4*(sqrt(2)*e^(-1/2*x^2)/sqrt(pi) - 1.41421356237309*e^(-0.500000000000000*x^2)/sqrt(pi))^2

which coincides with the answer of unpatched Sage.

Really, in some sense the new version is 'better' because it keeps the minus where the user put it. Do (-I - 1/3*sqrt(2))^2 and (I + 1/3*sqrt(2))^2 do anything different from before?

They do.

sage: (-I - 1/3*sqrt(2))^2
1/9*(-sqrt(2) - 3*I)^2
# Was: 1/9*(sqrt(2) + 3*I)^2
sage: (I + 1/3*sqrt(2))^2
1/9*(sqrt(2) + 3*I)^2
# Was: 1/9*(sqrt(2) + 3*I)^2

Really I don't understand where that comes from. Namely, as far as I know, the sign of gcd or lcm did not change. The only difference in gcd/lcm examples involving I and sqrt(2) is the fact that without my patch an error would be raised.

comment:25 in reply to: ↑ 24 Changed 9 years ago by kcrisman

Really, in some sense the new version is 'better' because it keeps the minus where the user put it. Do (-I - 1/3*sqrt(2))^2 and (I + 1/3*sqrt(2))^2 do anything different from before?

They do. sage: (-I - 1/3*sqrt(2))2 1/9*(-sqrt(2) - 3*I)2 # Was: 1/9*(sqrt(2) + 3*I)2

Okay, that's what I figured would happen.

sage: (I + 1/3*sqrt(2))2 1/9*(sqrt(2) + 3*I)2 # Was: 1/9*(sqrt(2) + 3*I)2

Okay, that is the same as before in any case.

Really I don't understand where that comes from. Namely, as far as I know, the sign of gcd or lcm did not change. The only difference in gcd/lcm examples involving I and sqrt(2) is the fact that without my patch an error would be raised.

Yes,

sage: gcd(I + 1/3*sqrt(2),I + 1/3*sqrt(2))
---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)

Interesting. Like I said earlier, Sage now sends such symbolic things to Pynac/Ginac?, which however needs to use Sage again in order to do certain calculations like this (Mike Hansen explained it well in your thread on sage-devel). So apparently in the past it decided to factor out a -1 in certain cases based on some gcd - maybe just a default of 1 or something - whereas now it does not, due to the new gcd. What does

sage: gcd(I + 1/3*sqrt(2),I + 1/3*sqrt(2))

do after the patch? What happens with this?

sage: gcd(4,2+2*I*sqrt(3))

Notice that I am purposely taking two things NOT in a polynomial or extension thing, where it can be calculated whether it's a UFD. In the symbolic ring, who knows what's "right"?

comment:26 Changed 9 years ago by dsm

I think I've figured out hg queues now, so I can actually try out the patch! [No guarantees, though.] Was fuzz-testing it and noticed the following:

sage: var("x")
x
sage: a = 1/3*x**0
sage: b = 2/5*x**0
sage: type(a)
<type 'sage.symbolic.expression.Expression'>
sage: gcd(a,b)
1/15
sage: lcm(a,b)
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)

/Applications/sage_devel/devel/sage-main/<ipython console> in <module>()

/Applications/sage_devel/local/lib/python2.6/site-packages/sage/rings/arith.pyc in lcm(a, b)
   1562                 return ZZ(a).lcm(ZZ(b))
   1563             except TypeError:
-> 1564                 raise TypeError, "unable to find lcm of %s and %s"%(a,b)
   1565         return LCM(b)
   1566 

TypeError: unable to find lcm of 1/3 and 2/5

All bets are off in SR, but there is an Expression.gcd (no .lcm, strangely). Now that QQ reduces to ZZ, maybe switch the attempt to try ZZ to QQ, so that rationals Just Work(tm)?

This surprised me a little too:

sage: R.<x>=QQ[]
sage: 
sage: a = 1/3*x**0
sage: b = 2/5*x**0
sage: type(a)
<type 'sage.rings.polynomial.polynomial_rational_flint.Polynomial_rational_flint'>
sage: gcd(a,b), lcm(a,b)
(1, 1)
sage: 
sage: a = (1*x**0)/(3*x**0)
sage: b = (2*x**0)/(5*x**0)
sage: type(a)
<class 'sage.rings.fraction_field_element.FractionFieldElement_1poly_field'>
sage: gcd(a,b), lcm(a,b)
(1, 1)

given that the seemingly harder cases

sage: a = 1/3+x
sage: b = 2/5+x
sage: type(a)
<type 'sage.rings.polynomial.polynomial_rational_flint.Polynomial_rational_flint'>
sage: gcd(a,b), lcm(a,b), a*b-gcd(a,b)*lcm(a,b)
(1, x^2 + 11/15*x + 2/15, 0)
sage: 
sage: a = (1*x**0)/(3*x**0)+x
sage: b = (2*x**0)/(5*x**0)+x
sage: type(a)
<class 'sage.rings.fraction_field_element.FractionFieldElement_1poly_field'>
sage: gcd(a,b), lcm(a,b), a*b-gcd(a,b)*lcm(a,b)
(1, x^2 + 11/15*x + 2/15, 0)

behave a little nicer. Is it possible to recover the rational results?

comment:27 follow-up: Changed 9 years ago by dsm

Ah, I just read the comments in #8111 (re: having to dive into flint). Well, at least changing ZZ to QQ in rings.arith.lcm would be straightforward.

comment:28 in reply to: ↑ 27 Changed 9 years ago by SimonKing

Replying to dsm:

Ah, I just read the comments in #8111 (re: having to dive into flint). Well, at least changing ZZ to QQ in rings.arith.lcm would be straightforward.

It could be that it is not straight forward. If my analysis of the segfault occuring with a previous version of my patch is correct, changing ZZ to QQ could very well trigger the segfault again.

comment:29 follow-up: Changed 8 years ago by dsm

  • Owner changed from AlexGhitza to dsm

Well, even after changing ZZ to QQ there everything seems to work as expected for me (no segfaults, a quick subset of tests all passed; I'll try testall long later and look for more crashers), and if it did explode it sort of feels like it'd be due to a bug somewhere else that it revealed.. but that can wait for another ticket later.

The real meat is already in your patch, and I'm already using it, so I'm very grateful for the work, and hope someone gives the okay soon! :^)

comment:30 Changed 8 years ago by dsm

  • Owner changed from dsm to AlexGhitza

comment:31 in reply to: ↑ 29 ; follow-up: Changed 8 years ago by SimonKing

Replying to dsm:

Well, even after changing ZZ to QQ there everything seems to work as expected for me (no segfaults, a quick subset of tests all passed; I'll try testall long later and look for more crashers), and if it did explode it sort of feels like it'd be due to a bug somewhere else that it revealed.. but that can wait for another ticket later.

What came out of the tests?

The real meat is already in your patch, and I'm already using it, so I'm very grateful for the work, and hope someone gives the okay soon! :^)

You (after finishing the tests, of course)?? Or you could post a patch for changing ZZ to QQ in rings.arith.lcm, and perhaps I can have a look on it? If I remember correctly, it is ok if you are reviewer for my patch and I am reviewer for your patch.

comment:32 in reply to: ↑ 31 Changed 8 years ago by dsm

Replying to SimonKing:

What came out of the tests?

I seem to recall it passed just fine but I can't find the output now :-( so I'm rerunning them (takes me about ~3h).

You (after finishing the tests, of course)??

Following was's recommendations (http://ask.sagemath.org/question/31/should-i-_really_-review-trac-tickets) I don't think I've got enough experience with the Sage codebase yet to set something to positive review; I'm still finding my way around.

comment:33 Changed 8 years ago by mariah

  • Reviewers changed from Marco Streng to Marco Streng, Mariah Lenox
  • Status changed from needs_review to positive_review

Applied patch to sage-4.7.1.alpha2, did 'sage -b' followed by 'make testlong'. All tests passed. Positive review!

comment:34 Changed 8 years ago by dsm

<Mr Burns> Excellent. </Burns> Many thanks to Simon for the work!

comment:35 Changed 8 years ago by jdemeyer

  • Milestone changed from sage-4.7.1 to sage-4.7.2
  • Type changed from defect to enhancement

comment:36 follow-up: Changed 8 years ago by jdemeyer

  • Status changed from positive_review to needs_work

Could you please change the commit message such that the lines are not so long? You should wrap long lines but make sure the first still makes sense by itself, that will appear in hg log.

Changed 8 years ago by SimonKing

Implement gcd and lcm for general fields, with the special case of the fraction field of a unique factorization domain

comment:37 in reply to: ↑ 36 Changed 8 years ago by SimonKing

  • Status changed from needs_work to positive_review

Replying to jdemeyer:

Could you please change the commit message such that the lines are not so long? You should wrap long lines but make sure the first still makes sense by itself, that will appear in hg log.

Admittedly the log line was "a bit" two long. It is now spread over three lines. I hope it is fine now.

comment:38 Changed 8 years ago by jdemeyer

  • Merged in set to sage-4.7.2.alpha0
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.