Opened 5 years ago
Closed 3 years ago
#15790 closed defect (fixed)
GCD of univariate polynomials: generic implementation for UFD, and sparse case
Reported by:  bruno  Owned by:  

Priority:  minor  Milestone:  sage6.9 
Component:  basic arithmetic  Keywords:  polynomial, gcd 
Cc:  vdelecroix  Merged in:  
Authors:  Bruno Grenet  Reviewers:  Julian Rueth 
Report Upstream:  N/A  Work issues:  
Branch:  d77ee90 (Commits)  Commit:  d77ee902e808082dcc44daf71bb0ea2756c621bb 
Dependencies:  Stopgaps: 
Description (last modified by )
[Update] Due to #13442, the purpose of this ticket changed slightly. It provides two functionalities:
 A generic implementation for the GCD of two polynomials over a UFD;
 A
gcd
method for sparse polynomials which allows to choose between thegeneric
implementation, or using the implementation fordense
polynomials. Default is generallygeneric
, but for polynomial overZZ
it isdense
. This is due to the very fast implementations available in Flint and NTL.
These implementations solve the problem mentioned in the first version of this ticket, that remains below for information.
[Old description, not relevant anymore]
The GCD of two sparse polynomials over ZZ
is not computed. If the polynomials are in dense representation or if they are in sparse representation but defined over QQ
, there is no problem.
sage: R1.<x>=PolynomialRing(ZZ,'x',sparse=true) sage: R2.<x>=PolynomialRing(ZZ,'x',sparse=false) sage: R3.<x>=PolynomialRing(QQ,'x',sparse=true) sage: p=x^23*x+2 sage: q=x^21 sage: gcd(R2(p),R2(q)) x  1 sage: gcd(R3(p),R3(q)) x  1 sage: gcd(R1(p),R1(q))  TypeError Traceback (most recent call last) <ipythoninput60d7380207fef> in <module>() > 1 gcd(R1(p),R1(q)) /home/bgrenet/bin/sage6.0x86_64Linux/local/lib/python2.7/sitepackages/sage/rings/arith.pyc in gcd(a, b, **kwargs) 1605 return ZZ(a).gcd(ZZ(b)) 1606 except TypeError: > 1607 raise TypeError("unable to find gcd") 1608 return GCD(b, **kwargs) 1609 TypeError: unable to find gcd
Change History (29)
comment:1 Changed 5 years ago by
comment:2 Changed 5 years ago by
I'd implement all 3 with an optional algorithm
argument since each has different strengths/weaknesses. For example, 2 is good as long as the polynomials are relatively dense, but for something like x^100000  2
, I wouldn't want to convert that to a dense polynomial because of the memory usage (and take the time to construct it). From that, I think 3 should be the default. Possible names:
"fraction_field"
"dense"
"pseudodivision"
Nevertheless, I think we should have some gcd
since this would be a surprise for the more casual user.
comment:3 Changed 5 years ago by
 Milestone changed from sage6.2 to sage6.3
comment:4 Changed 5 years ago by
 Cc vdelecroix added
comment:5 Changed 5 years ago by
 Branch set to u/bruno/gcd_of_sparse_univariate_polynomials_over_zz
comment:6 Changed 5 years ago by
 Commit set to 8ddbb9efad8d53552fa26727fb484720bda66025
Branch pushed to git repo; I updated commit sha1. New commits:
8ddbb9e  Add gcd method to class Polynomial_generic_sparse

comment:7 Changed 5 years ago by
I've implemented the gcd
method in the class Polynomial
using the pseudodivision algorithm, as well as a method gcd
in Polynomial_generic_sparse
which allows to choose the algorithm between "pseudodivision"
, "dense"
or "fraction_field"
.
Needs review!
One remaining limitation is the absence of quo_rem
method in Polynomial_generic_sparse
. I created a specific ticket #16544 for this.
comment:8 Changed 5 years ago by
 Keywords polynomial gcd added; GCD removed
 Status changed from new to needs_review
comment:9 Changed 5 years ago by
comment:10 Changed 5 years ago by
 Commit changed from 8ddbb9efad8d53552fa26727fb484720bda66025 to 65943e1a4ae6da9d2c7579ba14317aeada6eb102
Branch pushed to git repo; I updated commit sha1. New commits:
96f11a1  Add special case for p.pseudo_quo_rem(q) when q is a constant

9000fb8  Add quo_rem method in class Polynomial_generic_sparse

84833c9  Move inversion of unit

fb1d21c  Add doc tests

65943e1  Merge branch 't/16544/fb1d21cbee05d6be577bf66120aa20d3bd498ace' into t/15790/gcd_of_sparse_univariate_polynomials_over_zz

comment:11 Changed 5 years ago by
 Status changed from needs_review to needs_work
comment:12 Changed 5 years ago by
 Milestone changed from sage6.3 to sage6.4
comment:13 Changed 4 years ago by
 Commit changed from 65943e1a4ae6da9d2c7579ba14317aeada6eb102 to 5c84d7115a1e5406487b5d59746554b44be995e4
comment:14 Changed 4 years ago by
 Status changed from needs_work to needs_review
The ticket is ready for review. Note that as default, I have the "dense" algorithm for polynomials over ZZ and "pseudodivision" for the other ones. The idea is that arithmetic on ZZ[X] is very efficient.
Of course, gcd(x^1000001,x^10001)
is much faster (by a factor approx. 400) with "pseudodivision" algorithm because of memory allocation for the "dense" algorithm, and thanks to a very good behavior of the "pseudodivision" algorithm on this example. Yet, if you consider gcd(x^1000001,x^21)
(which should be easy btw), the "dense" algorithm is faster and I even stopped the "pseudodivision" algorithm after a quite long time.
Note: It may be possible to reduce the huge allocation time since with NTL implementation, this is much faster: https://groups.google.com/forum/#!topic/sagedevel/6qhW90dgd1k.
comment:15 followup: ↓ 17 Changed 4 years ago by
Do you have a heuristic for determining when one algorithm is faster than another (when the base ring is ZZ
)? It doesn't have to be a particularly good one, but considering the time differences you're getting, it seems like we should have some automatic choice of algorithm.
comment:16 Changed 4 years ago by
 Commit changed from 5c84d7115a1e5406487b5d59746554b44be995e4 to ec2ef11336040e95c8cdf1bed7523006747b5ab1
Branch pushed to git repo; I updated commit sha1. New commits:
ec2ef11  Choice between FLINT and NTL

comment:17 in reply to: ↑ 15 Changed 4 years ago by
Replying to tscrim:
Do you have a heuristic for determining when one algorithm is faster than another (when the base ring is
ZZ
)? It doesn't have to be a particularly good one, but considering the time differences you're getting, it seems like we should have some automatic choice of algorithm.
I don't have one. Based on the comments on https://groups.google.com/forum/#!topic/sagedevel/6qhW90dgd1k, I add the possibility¹ for the "dense" algorithm to use NTL instead of FLINT. After some tests, I chose a heuristic for using NTL rather than FLINT (basically, FLINT is faster for lowdegree polynomials and nottoosparse ones). Then the "dense" algorithm is almost always faster than "pseudodivision" but for some very specific examples (I haven't been able to build other ones): gcd(x^N1,x^M1)
with both N
and M
pretty large. An the ratio is not hundreds anymore, but approx. 10.
¹ Actually, no choice for the user but only an automatic switch.
comment:18 Changed 4 years ago by
 Status changed from needs_review to needs_work
needs work, does not build, see patchbot report:
File "sage/structure/element.pyx", line 3433, in sage.structure.element.NamedBinopMethod.__init__ (build/cythonized/sage/structure/element.c:26905) AttributeError: 'sage.structure.element.NamedBinopMethod' object has no attribute '__name__'
comment:19 Changed 4 years ago by
 Commit changed from ec2ef11336040e95c8cdf1bed7523006747b5ab1 to 3015456844cb54f37cf4bbae8423521cdd33624a
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
b272960  Add pseudodivision and gcd for Polynomials (as generic as possible)

c88b600  Add gcd method to class Polynomial_generic_sparse

6f6a167  Add special case for p.pseudo_quo_rem(q) when q is a constant

9f797ef  Correct formatting of doctests + raise error for unknown algorithm

6c31b67  Change Error for pseudodivision by zero

f623b88  Choice between FLINT and NTL

3015456  Adapted to 6.9beta1 => Merged two gcd methods

comment:20 Changed 4 years ago by
 Description modified (diff)
 Milestone changed from sage6.4 to sage6.9
 Status changed from needs_work to needs_review
 Summary changed from GCD of sparse univariate polynomials over ZZ to GCD of [sparse] univariate polynomials [over ZZ] [updated]
This is rebased, and now works (as far as I can tell).
comment:21 Changed 4 years ago by
 Status changed from needs_review to needs_work
Doc does not build
OSError: [polynomia] docstring of sage.rings.polynomial.polynomial_element.Polynomial.pseudo_quo_rem:38: WARNING: Duplicate explicit target name: "gtm138".
comment:22 Changed 3 years ago by
 Commit changed from 3015456844cb54f37cf4bbae8423521cdd33624a to 47badc1b7322077ffbf2f89ab6ad2acc2dce18b4
Branch pushed to git repo; I updated commit sha1. New commits:
47badc1  Move fallback algorithm to category UFD

comment:23 Changed 3 years ago by
 Commit changed from 47badc1b7322077ffbf2f89ab6ad2acc2dce18b4 to f9089369482cc589372b0f3a7461d633d9540964
Branch pushed to git repo; I updated commit sha1. New commits:
f908936  Simplify GCD for sparse polynomials

comment:24 Changed 3 years ago by
 Description modified (diff)
 Keywords sparse removed
 Status changed from needs_work to needs_review
 Summary changed from GCD of [sparse] univariate polynomials [over ZZ] [updated] to GCD of univariate polynomials: generic implementation for UFD, and sparse case
In the two latest commits, I propose the following changes:
 Instead of having a fallback algorithm in the method
Polynomial.gcd
which tests whether the base ring in a unique factorization domain, this algorithm is moved to the categoryUniqueFactorizationDomains
.  I simplified the method
Polynomial_generic_sparse.gcd
: There remain only two possibilities (out of three), namely the one previously calledpseudodivision
(renamedgeneric
, since it uses the generic code) and thedense
version, which may in some cases be faster, for instance overZZ
. I removed thefraction_field
option which does not make much sense I think.
At the same time, I also take Frédéric's comment into account.
comment:25 Changed 3 years ago by
 Branch changed from u/bruno/gcd_of_sparse_univariate_polynomials_over_zz to u/saraedum/gcd_of_sparse_univariate_polynomials_over_zz
comment:26 Changed 3 years ago by
 Commit changed from f9089369482cc589372b0f3a7461d633d9540964 to d77ee902e808082dcc44daf71bb0ea2756c621bb
Branch pushed to git repo; I updated commit sha1. New commits:
d77ee90  Improved docstring of PolynomialElement.pseudo_quo_rem()

comment:27 Changed 3 years ago by
 Status changed from needs_review to positive_review
Made some minor docstring changes. Everything else looks fine.
comment:28 Changed 3 years ago by
 Reviewers set to Julian Rueth
Thanks for improving the docstring, and for the review! I add your name as a reviewer.
comment:29 Changed 3 years ago by
 Branch changed from u/saraedum/gcd_of_sparse_univariate_polynomials_over_zz to d77ee902e808082dcc44daf71bb0ea2756c621bb
 Resolution set to fixed
 Status changed from positive_review to closed
For
PolynomialRing(ZZ,'x',sparse=false)
, one of FLINT or NTL is used (FLINT by default).For
PolynomialRing(QQ,'x',sparse=true)
, thegcd
computation is based onquo_rem
. This is not possible over the integer sincequo_rem
works for polynomials over a field.There are several possible solutions:
PolynomialRing(QQ,'x',sparse=true)
; compute their GCDg
, the LCMl
of the denominators of the coefficients ofg
, and returnl*g
.PolynomialRing(ZZ,'x',sparse=false)
; compute their GCD and return the result in sparse representation with a second coercion.Solution 3 seems cleaner but does not give especially efficient results since the algorithms described in Cohen GTM 138 have complexity polynomial in the degree of the input polynomial. Actually, computing the GCD of two sparse polynomial is NPhard ![1] so no better solution than converting to the dense representation is likely to exist. To my mind, Solution 2 is then the best one.
Yet, I am new to Sage so comments are welcome!
![1] D. Plaisted, Sparse complex polynomials and polynomial reducibility. J. Comput. Syst. Sci. 14 (2), 210–221, 1977.