fmpz_mpoly.gcd does not include content of polynomials

[RotemKalisch]

When using fmpz_mpoly, the gcd contains the content of the polynomials:

>>> import flint
>>> ctx = flint.fmpq_mpoly_ctx.get('n')
>>> flint.fmpz_mpoly(4, ctx).gcd(flint.fmpz_mpoly(6, ctx)) 
2

However, when using fmpq_mpoly, gcd does not contain the content of the polynomials:

>>> import flint
>>> ctx = flint.fmpq_mpoly_ctx.get('n')
>>> flint.fmpq_mpoly(4, ctx).gcd(flint.fmpq_mpoly(6, ctx)) 
1

Digging further, I've noticed that fmpq does not have a gcd method, which could explain this root cause.
The GCD of two rationals a/b, c/d is simply extended as gcd(a*d, b*c) / (b*d).

I'm not sure if this falls under the category of a bugfix or a feature request, but is it possible to achieve this behavior?

[oscarbenjamin]

The mpoly types have not been included in a full release yet (I guess you are using a prerelease?) so it should still be possible to change this without breaking compatibility.

Currently this just calls the FLINT function fmpq_mpoly_gcd and returns the result as is. Looks like it always makes the GCD monic:

In [6]: [n] = ctx.gens()

In [7]: (4*(n**2-1)).gcd(2*(n-1))
Out[7]: n - 1

It would be worth reviewing all poly/mpoly types to see if this is consistent.

[oscarbenjamin]

Related #189

[RotemKalisch]

I'm indeed using a pre-release. Specifically version 0.7.0a5.

It looks like the decision to return only the monic GCD was made in FLINT.
However, FLINT also provides a method to extract the content of the GCD - fmpq_mpoly_content, which is not present in the fmpq_mpoly.

To me the decision to exclude the content from the GCD is counter-intuitive, but if you wish to maintain compatibility with FLINT, wrapping the fmpq_mpoly.content() as well as wrapping a fmpq.gcd() method should allow users to achieve a full GCD functionality.

[oscarbenjamin]

It seems to be consistent with fmpq_poly and nmod_poly and the corresponding mpoly types e.g.:

In [2]: x = flint.fmpq_poly([0, 1])

In [3]: (4*(x**2-1)).gcd(2*(x-1))
Out[3]: x + (-1)

In [4]: x = flint.nmod_poly([0, 1], 11)

In [5]: (4*(x**2-1)).gcd(2*(x-1))
Out[5]: x + 10

In [27]: ctx = flint.fmpz_mod_mpoly_ctx.get('n', 11)

In [28]: [n] = ctx.gens()

In [29]: (-4*(n**2-1)).gcd(-2*(n-1))
Out[29]: n + 10

I guess it is natural for GCD of a polynomial over a field to be monic.

[oscarbenjamin]

It would be reasonable top add a content method. There is term_content but that is the GCD of the monomials rather than the coefficients:

In [33]: (4*n**3 + 2*n**2).term_content()
Out[33]: n^2

[oscarbenjamin]

Actually term content is the GCD of the terms:

In [42]: ctx = flint.fmpz_mpoly_ctx.get('n')

In [43]: [n] = ctx.gens()

In [44]: (4*n**3 + 2*n**2).term_content()
Out[44]: 2*n^2

Just when the ground domain is a field the gcd of the terms is defined as always monic.

[oscarbenjamin]

In general it would be good to have primitive and content for all polynomials. Currently only fmpz_mpoly has these. It doesn't look like gr_poly and gr_mpoly have anything like this.

Returning a monic GCD seems to be consistent with SymPy's analogous types over field domains:

In [7]: R = QQ[x]

In [8]: xR = R(x)

In [9]: R.gcd(4*(xR**2-1), 2*(xR-1))
Out[9]: x - 1

In [10]: R = ZZ[x]

In [11]: xR = R(x)

In [12]: R.gcd(4*(xR**2-1), 2*(xR-1))
Out[12]: 2*x - 2

In [13]: Poly(4*(x**2-1), x, domain=QQ).gcd(Poly(2*(x-1), x, domain=QQ))
Out[13]: Poly(x - 1, x, domain='QQ')

In [14]: Poly(4*(x**2-1), x, domain=ZZ).gcd(Poly(2*(x-1), x, domain=ZZ))
Out[14]: Poly(2*x - 2, x, domain='ZZ')

In [25]: R = QQ.frac_field(y)[x]

In [26]: R
Out[26]: QQ(y)[x]

In [27]: R.gcd((y**2*xR**2-1), (y*xR-1))
Out[27]: x - 1/y

[oscarbenjamin]

I'm not sure that content is what you want here. What is removed from the GCD is the leading coefficient. What you want I think is some GCD of the leading coefficients:

In [8]: p1 = 2*(n**2 - 1)

In [9]: p2 = 2*(n - 1)

In [10]: p1.gcd(p2)
Out[10]: n - 1

In [11]: p1.leading_coefficient()
Out[11]: 2

In [12]: p2.leading_coefficient()
Out[12]: 2

In [13]: c1 = p2.leading_coefficient()

In [14]: c2 = p2.leading_coefficient()

In [15]: c1.numer().gcd(c2.numer())
Out[15]: 2

[RotemKalisch]

Seems like you are correct.

This example displays this behavior even better:

In [36]: p1 = flint.fmpq_mpoly(4*n**2 + 2*n, ctx)

In [37]: p2 = flint.fmpq_mpoly(2*n**2 + n, ctx)

In [38]: p1.gcd(p2)
Out[38]: n^2 + 1/2*n

In [39]: p1.leading_coefficient()
Out[39]: 4

In [40]: p1.leading_coefficient().numer().gcd(p2.leading_coefficient().numer())
Out[40]: 2

Still, one can run

In [41]: p3 = flint.fmpq_mpoly(n**2 / 2 + n / 4, ctx)

In [42]: p1.gcd(p3)
Out[42]: n^2 + 1/2*n

And in this case we cannot simply use p3.leading_coefficient().numer(), because we're dealing with a rational number.
So I'd still argue that fmpq.gcd(fmpq) is needed.

Either way - the solution using the numer() part solves my specific problem. Thank you for your help :)

[oscarbenjamin]

So I'd still argue that fmpq.gcd(fmpq) is needed.

You can still implement it as you suggested before.

This one is provided in FLINT though as fmpq_gcd and is also consistent with SymPy's behaviour:

In [1]: QQ.gcd(QQ(2,3), QQ(4,9))
Out[1]: 2/9

The doc for fmpq_gcd notes that it is inconsistent with fmpq_poly_gcd.

I'm not sure though if it would be better to change the definition of gcd for fmpq_poly and fmpq_mpoly. The generics use monic GCD it seems:

In [7]: from flint.types._gr import gr_fmpq_ctx, gr_gr_poly_ctx

In [8]: ctx = gr_gr_poly_ctx.new(gr_fmpq_ctx)

In [9]: x = ctx.gen()

In [10]: (2*x).gcd(4*x)
Out[10]: x

[oscarbenjamin]

I've added fmpq.gcd in gh-250.