Opened 12 years ago

Closed 11 years ago

# gcd and lcm for fraction fields

Reported by: Owned by: Simon King Alex Ghitza major sage-4.7.2 basic arithmetic gcd lcm fraction fields Burcin Erocal sage-4.7.2.alpha0 Simon King Marco Streng, Mariah Lenox N/A

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.

### comment:1 Changed 12 years ago by Simon King

Description: modified (diff)

### comment:2 Changed 12 years ago by Simon King

Authors: → Simon King new → needs_review PLEASE CHANGE → defect

```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:  6 Changed 12 years ago by Marco Streng

Reviewers: → Marco Streng

• 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 12 years ago by Luis Felipe Tabera Alonso

Related tickets: #9819, #10459

### comment:5 Changed 12 years ago by D.S. McNeil

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

### comment:6 in reply to:  3 Changed 12 years ago by Simon King

Hi Marco,

• 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:  11 Changed 12 years ago by Simon King

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:  9 Changed 12 years ago by Simon King

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 12 years ago by Simon King

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:  22 Changed 12 years ago by Marco Streng

Status: needs_review → 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:  12 Changed 12 years ago by Luis Felipe Tabera Alonso

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?

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 12 years ago by Simon King

Hi Luis!

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:  14 Changed 12 years ago by Marco Streng

• 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 12 years ago by Simon King

• 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 12 years ago by Simon King

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 12 years ago by Simon King

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:  18 Changed 12 years ago by Simon King

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 12 years ago by Simon King

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 12 years ago by Simon King

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 12 years ago by Karl-Dieter Crisman

Just FYI, #8111 is possibly also related.

### comment:21 follow-up:  23 Changed 12 years ago by Simon King

Cc: Burcin Erocal added needs_work → 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 12 years ago by Simon King

Hi Marco,

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:  24 Changed 12 years ago by Karl-Dieter Crisman

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:  25 Changed 12 years ago by Simon King

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 12 years ago by Karl-Dieter Crisman

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 12 years ago by D.S. McNeil

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:  28 Changed 12 years ago by D.S. McNeil

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 12 years ago by Simon King

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:  31 Changed 12 years ago by D.S. McNeil

Owner: changed from Alex Ghitza to D.S. McNeil

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 12 years ago by D.S. McNeil

Owner: changed from D.S. McNeil to Alex Ghitza

### comment:31 in reply to:  29 ; follow-up:  32 Changed 12 years ago by Simon King

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 12 years ago by D.S. McNeil

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 11 years ago by Mariah Lennox

Reviewers: Marco Streng → Marco Streng, Mariah Lenox needs_review → 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 11 years ago by D.S. McNeil

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

### comment:35 Changed 11 years ago by Jeroen Demeyer

Milestone: sage-4.7.1 → sage-4.7.2 defect → enhancement

### comment:36 follow-up:  37 Changed 11 years ago by Jeroen Demeyer

Status: positive_review → 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 11 years ago by Simon King

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 11 years ago by Simon King

Status: needs_work → positive_review

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`.